
Calhoun 

iniutuiiaiul AKliiv« ou tfit Nilvdl Poi($ra{jua(« School 


Calhoun: The NPS Institutional Archive 
DSpace Repository 



Theses and Dissertations 


1. Thesis and Dissertation Collection, all items 


1994-12 

Kinematics and dynamics of a cyclonic eddy 
off Pt. Arena, California 

Chumbinho, Rogerio Paulo Antunes 

Monterey, California. Naval Postgraduate School 


http://hdl.handle.net/10945/42791 


This publication is a work of the U.S. Government as defined in Title 17, United 
States Code, Section 101. Copyright protection is not available for this work in the 
United States. 

Downloaded from NPS Archive: Calhoun 



DUDLEY 

KNOX 

LIBRARY 


http://www.nps.edu/ljbrary 


CsMwun is the Neval Postgraduate School's public access distal repository for 
research oiateriels and tnstitutjiooal pubftcatiions created by the NPS community. 
Cathouni is named for Professor of Mathematcs Guy K. CatHiuo, NPS's first 
appointed — and publi^d — scholar^ author. 

Dudley Knox Library / Naval Postgraduate School 
411 Dyer Road / 1 University Circle 
MontereVr California USA 93943 



NAVAL POSTGRADUATE SCHOOL 
Monterey, California 



^Original oontalns eibloi»'s^ 
plates: All DSIC veproduot* 
Ions Till be in bi a o fc - 
ThlU*-' ■ 


DISSERTATION 

KINEMATICS AND DYNAMICS OF A CYCLONIC 
EDDY OFF PT. ARENA, CALIFORNIA 

by 

Rogerio Paulo Antunes Chumbinho 
December 1994 

Thesis Advisor: Robert L. Haney 

Approved for public release; distribution is unlimited 


19950830 055 


DTI® QOALrr? 








REPORT DOCUMENTATION PAGE 

Form Approved 

0MB No. 0704-0188 

Public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, 
gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this 
collection of information, including suggestions for reducing this burden, to Washington Headquarters Services, Directorate for Information Operations and Reports, 1215 Jefferson 
Davis Highway, Suite 1204, Arlington, VA 22202^302, and to the Office of Management and Budget, Paperwork Reduction Project (0704-0188), Washington, DC 20503. 

1. AGENCY USE ONLY (Leave blank) 2. REPORT DATE 

December 1994 

3. REPORT TYPE AND DATES COVERED 

Doctoral Dissertation 

4. TITLE AND SUBTITLE 

KINEMATICS AND DYNAMICS OF A CYCLONIC EDDY OFF Pt. ARENA, 
CALIFORNIA 

5. FUNDING NUMBERS 

6. AUTHOR(S) 

Rog6rio Paulo Antunes Chumbinho 

7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 

Naval Postgraduate School 

Monterey, CA 93943-5000 

8. PERFORMING ORGANIZATION 

REPORT NUMBER 

9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESSE(S) 

10. SPONSORING/MONITORING 

AGENCY REPORT NUMBER 

11 . SUPLEMENTARY NOTES The views expressed in this thesis are those of the author and do not reflect the official policy 
or position of the Department of Defense or the U.S. Government. 

12a. DISTRIBUTION/AVAILABILITY STATEMENT 

Approved for public release; distribution is unlimited 

12b. DISTRIBUTION CODE 

13. ABSTRACT (maximum 200 words) 

In situ hydrographic data and remotely sensed sea surface temperature data were analyzed to describe the 
kinematics and dynamics of a cyclonic eddy observed off Point Arena, California, in May 1993. The hydrographic 
data were first objectively an^yzed using multiquadric-biharmonic interpolation (MQ-B) to obtain a three- 
dimensional field of density. This was then used in a primitive equation (PE) model with a digital filter initialization 
method to arrive at the dynamically balanced three-dimensional velocity field, including the vertical motion. Several 
aspects related to geostrophic adjustment, to data analysis and to the lack of exact synopticity in the observations 
were examined; these include the effects of internal tides, eddy rotation and MQ-B parameters sensitivity analysis. 

The results reveal a cyclonic circulation with meanders in the otherwise nearly circular flow, and maximum 
horizontal speeds of about 40-50 cm/s at the surface. The meanders have associated patterns of radially aligned 
patches of sinking and rising motion, extending coherently to about 500 m depth, with maximum vertical velocities 
of 20-22 m/day between 100 and 150 m depth. Because of the strong horizontal currents and the short horizontal 
length scales involved, the particle vertical excursions are very small. 

14. SUBJECT TERMS 

Vertical Velocity, Primitive Equation Modeling, Digital Filter Initialization, Quasigeo- 
strophic Theory, Omega Equation, Multiquadric-Biharmonic Interpolation. 

15. NUMBER OF PAGES 

98 

16. PRICE CODE 

17. SECURITY CLASSIFICATION 18. SECURITY CLASSIFICATION 

OF REPORT OF THIS PAGE 

Unclassified Unclassified 

19. SECURITY CLASSIFICATION 

OF ABSTRACT 

Unclassified 

19. LIMIITATION OF ABSTRACT 

UL 


NSN 7540-01-280-5500 


Standard Form 298 (Rev. 2-89) 
Prescribed by ANSI Std. Z39-18 
298-102 






















UNCLASSIFIED 

SECURITY CLASSIFICATION OF THIS PAGE 


13. ABSTRACT (continued) 

The vertical motion in the domain was also diagnosed from the analyzed density field using the quasigeo- 
strophic (QG) omega equation. The resulting vertical velocities are in good agreement with those from the initial¬ 
ized PE model except that the QG values were about 30 percent stronger in magnitude. These discrepancies were 
shown to be caused by the neglect of ageostrophic advection in the QG system. 


Accesion For 

NTIS CRA&I 
OTIC TAB 
Unannounced 
Justification 


□ 

□ 


By_ 

Distribution / 


Availability Cedes 



Avail and/or 

Dist 

Special 

/)'/ 



Standard Form 298, (Reverse) 


ii 


SECURITY CLASSIFICATION OF THIS PAGE 

UNCLASSIFIED 






Approved for public release; distribution is unlimited 


KINEMATICS AND DYNAMICS OF A CYCLONIC EDDY 
OFF PT. ARENA, CALIFORNIA 


Rog^rio Paulo Antunes Chumbinho 
Lieutenant, Portuguese Navy 
M.S., Naval Postgraduate School, 1993 

DOCTOR OF PHILOSOPHY IN PHYSICAL OCEANOGRAPHY 

from the 

NAVAL POSTGRADUATE SCHOOL 
December 1994 


Author: 
Approved by: 



Rog6rio Paulo Antunes Chumbinho 



Peter C. Chu 

Professor of Oceanography 


teven R. Ramp 
Professor of Oceanography 



Approved by: 


Richard S. Elster, Dean of Instruction 















ABSTRACT 


In situ hydrographic data and remotely sensed sea surface temperature data were 
analyzed to describe the kinematics and dynamics of a cyclonic eddy observed off Point 
Arena, California, in May 1993. The hydrographic data were first objectively analyzed 
using multiquadric-biharmonic interpolation (MQ-B) to obtain a three-dimensional field 
of density. This was then used in a primitive equation (PE) model with a digital filter 
initialization method to arrive at the djmamically balanced three-dimensional velocity 
field, including the vertical motion. Several aspects related to geostrophic adjustment, to 
data analysis and to the lack of exact synopticity in the observations were examined; these 
include the effects of internal tides, eddy rotation and MQ-B parameters sensitivity 
analysis. 

The results reveal a cyclonic circulation with meanders in the otherwise nearly circular 
flow, and maximum horizontal speeds of about 40-50 cm/s at the surface. The meanders 
have associated patterns of radially aligned patches of sinking and rising motion, extending 
coherently to about 500 m depth, with maximum vertical velocities of 20-22 m/day 
between 100 and 150 m depth. Because of the strong horizontal currents and the short 
horizontal length scales involved, the particle vertical excursions are very small. 

The vertical motion in the domain was also diagnosed from the analyzed density field 
using the quasigeostrophic (QG) omega equation. The resulting vertical velocities are in 
good agreement with those fi’om the initialized PE model except that the QG values were 
about 30 percent stronger in magnitude. These discrepancies were shown to be caused by 
the neglect of ageostrophic advection in the QG system. 






The views expressed in this dissertation are those of the author and do not 
reflect the official policy or position of the Department of Defense or the 
U. S. Government. 
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1. INTRODUCTION 


The large scale and the mesoscale circulation in eastern boundary currents (EBC) 
are becoming well established areas of active research. This is particularly true in the 
coastal zone, where the interaction between the eastern boundary currents and the boundary 
is strongest. Consequently, a series of research experiments designed to study the dynamics 
of these regions are being carried out. An area that has seen a great deal of research is off 
the west coast of the United States, mainly the Oregonian and Northern Californian coasts 
(California Cooperative Oceanic Fisheries Investigation, Coastal Ocean Dynamics 
Experiment, Coastal Transition Zone experiments, among others). 

An extensive work describing the circulation of the California Current System, 
assembling results from various early experiments, is found in Hickey (1979). This study 
documents the existence of a broad, slow southward flow offshore, with strong coastal 
up welling in selected areas during periods of strong northerly wind events. The study also 
describes an undercurrent flowing poleward and the reversal of the inshore surface current, 
that flows poleward when the wind relaxes. 

The early interest in eastern boundary currents was indeed linked to the study of 
upwelling, an important phenomenon that takes place over large areas of western 
continental margins and which contributes to the enhancement of the nutrient content of 
coastal waters (e.g. Huyer, 1983 and Huyer, 1984). It was later recognized, mainly with the 
availability of satellite imagery (Rienecker and Mooers, 1989; Pares-Sierra et at, 1993), 
that eastern boundary currents, once considered slow and steady flows, are in fact very rich 
in mesoscale features mainly resulting from the interaction of the currents with the coast, 
bottom, and local winds. From this interaction there results an intricate circulation where 
cold filaments extending offshore and permanent and recurring eddies are common 
features. 

With the identification of these mesoscale features, the ocean sciences community 
started to investigate reasons for their formation and development. This led to a wealth of 
observational programs utilizing a variety of remotely sensed and in situ techniques, not 
only off the west coast of the United States but also in other eastern boundary current 
regions (e.g.. Ramp et al, 1991; Sousa, 1992; Le Groupe Tourbillon, 1983). The 
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importance of mesoscale eddies in the transport of heat, salt and other bio-physical 
properties in the coastal ocean is also being investigated. 

It now seems clear that coastal upwelling, mesoscale eddies and cold filaments are 
closely related. An equatorward coastal jet will develop meanders by virtue of instabilities 
in the flow, associated with vorticity and curvature mechanisms. The mean relative 
vorticity of the jet increases toward the coast, and the ocean depth decreases in this 
direction, assuming a sloping bottom. This causes the mean potential vorticity to increase 
towards the coast, so the meandering of the jet behaves like a planetary Rossby wave, 
propagating upstream relative to the flow. This propagation speed can be of the same order 
of magnitude as the mean equatorward jet speed, and thus some meanders can become 
stationary. This effect, in conjunction with local topographic features such as capes, which 
tend to induce their own vorticity effects in the flow and also cause appreciable local wind 
changes, is a good scenario for the appearance of permanent or recurrent eddies at certain 
locations along the coast. These mechanisms may produce the offshore surface flow of 
recently upwelled waters, that tends to be stronger just south of the capes. In this situation, 
an eddy could entrain and advect offshore the upwelled coastal waters. 

The eddies will not remain stationary once generated. According to the analytical 
theory of Nof (1981,1984), the eddies will propagate westward, away from the boundary. 
Although Nof s analyses are based on somewhat strong assumptions, the general results are 
valid in a qualitative point of view and have been shown to be in good agreement with some 
physical cases (Bowman, 1985; Andrews and Scully-Power, 1976). This westward 
propagation could contribute to the formation of the filaments observed in the vicinity of 
capes; the filaments often have two or more coupled eddies at their extremity. 

A useful technique for the study of eddies and cold filaments is numerical 
simulation (e.g. Haidvogel et ah, 1991; Batteen etal, 1992; Pares-Sierra et al., 1993). In 
the California Current region, Paduan and Niiler (1990) and Swenson et al. (1992) applied 
quasigeostrophic theory to the motion of clusters of surface drifters. Their results suggest 
that meanders in the seaward flowing jets and filaments have secondary circulations and 
vertical motions of 10-20 m/day at the base of the surface mixed layer. Walstad et al. 
(1991) used a quasigeostrophic model to study the California coastal jet observed during 
the 1987 Coastal Transition Zone (CTZ) program. The focus of this study was to determine 
the rotational (nondivergent) part of the currents including the “level of no motion” which 
was found to be shallower toward the coast. It is important to note that all of the previous 
efforts to analyze the structure and dynamics of the California coastal jet and eddies have 
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been based on quasigeostrophic theory, and most of the studies such as those based on 
surface drifters, have been restricted to only the upper mixed layer. 

Other authors have used quasigeostrophic theory to analyze similar situations in 
other ocean areas. In the study of vertical motion near the Gulf Stream, Lindstrom and 
Watts (1994) used the quasigeostrophic vorticity equation, among other methods, to 
compute vertical motion making direct use of observed data, since their data had sufficient 
temporal resolution. Henrick et al. (1979) made use of an analytical model based on the 
quasigeostrophic potential vorticity balance to predict perturbations of pressure, density 
and velocities in mesoscale eddies. Arhan and Colin de Verdiere (1985) used data collected 
during the Tourbillon experiment to do a statistical fit of streamfunction and 
quasigeostrophic vorticity fields, and then to describe the dynamics of an anticyclonic 
eddy. The vertical motion in this case was not well described. A number of other studies 
present only a description of data collected in the vicinity of eddies, without attempting to 
go beyond basic analyses. 

The study described here deals with such an eddy observed offshore of Pt. Arena, 
California. The eddy is about 90 km in diameter and 1000 m depth. The eddy was observed 
between 22 and 25 May, 1993, during a cruise on board the R/V Point Sur. The processes 
described above are particularly significant in the area where the eddy was observed; the 
region between Pt. Arena and Pt. Reyes has the most favorable upwelling conditions of the 
entire coast (Huyer, 1984). A study of the California Current System off Oregon and 
Northern California by Dceda and Emery (1984) includes satellite imagery showing 
filaments and eddies offshore of Pt. Arena in the summer of 1982, suggesting the recurrent 
nature of the eddy observed in May 1993. 

The objective of this study is to describe the complete three-dimensional 
circulation associated with a closed cyclonic eddy that was observed in the California 
coastal region. The description is accomplished by assimilating a unique data set into a 
primitive equation numerical model. The methodology used is to objectively analyze the 
quasi-synoptic hydrographic data collected during a recent Eastern Boundary Current 
(EBC) research program, and satellite sea surface temperature data, and then to use the 
analyzed density field as input to the model. The use of a new initialization procedure yields 
dynamically balanced velocity and mass fields that are used in the description of the eddy 
and its dynamics. This methodology is not restrictive to the application above (description 
of the observed cyclonic eddy) or to the California coast. It can be applied in other areas of 
the world ocean where mesoscale vertical motion has been poorly or not determined at all. 


3 



as long as the available data has enough spatial resolution. These areas include, for 
example, other eastern boundary current regions (e.g., the Eastern Atlantic Current along 
the coast of Portugal). Other possible applications are the study of the three-dimensional 
circulation associated with coastal jets and filaments, and frontal zones. 

The interpretation of the results from the model includes a discussion of aliasing or 
distortion in the data caused by processes with temporal and spatial scales of the same order 
of magnitude of the spacing of the hydrographic stations. Other important questions 
addressed in this study are whether internal tides were present during the collection of the 
hydrographic data, and what their effect might be, and also the distortion effects related to 
the rotation of the cyclonic eddy during the CTD sampling. 

In addition to presenting a complete description of the eddy, the hypothesis that the 
dynamics of the eddy can be interpreted within the framework of quasigeostrophic theory 
was also tested; quasigeostrophic theory describes adequately the dynamics of systems 
with small Rossby numbers. Attention was focused on the vertical motion, for several 
reasons. First, it is an important quantity that has not been measured directly in the EEC 
program; secondly, it is responsible for the vertical advection of important tracers and it 
influences biological productivity; thirdly, it is related to departures from geostrophy and 
secondary circulations; finally, it is essential in some conceptual models of vorticity in 
meandering jets and eddies. 

This study is organized as follows: Chapter n has two parts; in the first part the 
available data is described, including hydrographic and moored data; the second part 
presents the data analysis, including the interpolation to a regular grid for model 
applications. Chapter El describes the methods used for the computation of the vertical 
velocity, namely the circulation model initialized in such a way as to isolate the slow mode 
captured by the data, and the quasigeostrophic omega equation. Chapter IV has two parts: 
the first part presents and discusses the results obtained with these methods and describes 
the physical features of the eddy. It also includes discussions on the initialization and 
geostrophic adjustment, and potential data related problems such as the effect of internal 
tides, eddy rotation and sensitivity to data interpolation. The second part presents and 
discusses the results obtained with the quasigeostrophic omega equation. Finally, 
Chapter V gives conclusions and final recommendations. 
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II. DATA 


A. DATA DESCRIPTION 

The data used in this study consists of satellite data and hydrographic data from a 
quasi-synoptic Eastern Boundary Current (EBC) CTD survey off Pt. Arena, California, in 
May 1993. The objective of the survey was to observe and resolve eddies and filaments that 
had been seen off Pt, Arena in earlier studies such as the Coastal Transition Zone (CTZ) 
program in 1987-1988. As shown below, the survey was successful in mapping a very large 
fraction of a closed cyclonic eddy that had just cut off and separated from the coast, 

1. Satellite imagery 

The satellite data used in this work consists of images obtained at the Naval 
Postgraduate School, The raw images contain all five chaimels of the Advanced Very High 
Resolution Radiometer (AVHRR) flying on the National Oceanic and Atmospheric 
Administration (NOAA) series of Earth observing satellites. This sensor observes the 
surface of the Earth in the visible and infrared (IR) regions, using wavelengths between 3.7 
and 11.5 pm. The sea surface temperatures (SST) were computed by means of a 
multispectral algorithm that converts brighmess temperatures to sea surface temperature 
and corrects for the effect of the atmospheric water vapor on the observed radiances. 

The images were processed using software developed at the University of Miami. 
This processing consisted of correcting the images for appropriate navigation 
(geographical correction) and of remapping the entire image into the region of interest, 
centered at a point off Point Arena, California. The ground resolution depends on the 
distance between each particular pixel and the satellite nadir and is estimated in the range 
of 1 to 3 km. 

Most of the images collected during the period of the cruise contained large 
amounts of cloud cover; for this reason, only the images for 21 and 22 May 1993 (23:00 
Greenwich Mean Time) were suitable to be used in this work. An example of an image 
from 22 May, obtained from only the channel four of the AVHRR, containing also 
information on the area covered by the model domain (see Chapter HI) and hydrographic 
stations (see section n.A.2), is given in Figure 1. 
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Figure 1. National Oceanic and Atmospheric Administration (NOAA) satellite image, 
obtained from Channel 4 of the Advanced Very High Resolution Radiometer 
(AVHRR), of the Sea Surface Temperature (SST) off Point Arena, on 22 May 1993, 
23 GMT. The crosses represent the positions of the CTD stations: also shown is the 
area of study (model domain). Temperatures are in degrees Celsius. 
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A surface temperature map was obtained by sampling the multispectral image of 
May 22, 23:00 GMT at the same horizontal positions of the grid used in the model, thus 
providing a grid of satellite-derived SST used later in the data interpolation (see 
section n.B.a). 

2. Hydrographic data 

The hydrographic data used in this study were collected during a Naval 
Postgraduate School cruise on board R/V Point Sur, between 21 and 25 May, 1993, using 
a Neil Brown Mark El Conductivity, Temperature and Depth (CTD) instrument. The CTD 
casts were performed at an array of stations listed in Table 1. The column labeled MAX Z 
means the maximum depth attained in the cast; the column labeled DEPTH is the depth to 
the bottom at the station. See also Figures 1 and 2. 


LONG (W) 

LAT(N) 

StaNo 

DATE 

TIME (h) 

MAXZ 

(m) 

DEPTH 

(m) 

-123.697 

38.546 

1 

22 

23.9 

479 

500 

-123.779 

38.643 

2 

22 

22.7 

419 

460 

-123.857 

38.750 

3 

22 

21.4 

371 

410 

-123.951 

38.859 

4 

22 

20.0 

415 

455 

-124.027 

38.970 

5 

22 

18.5 

471 

530 

-124.176 

38.889 

6 

23 

10.6 

2253 

2420 

-124.080 

38.782 

7 

23 

8.3 

2041 

2110 

-123.984 

38.675 

8 

23 

6.0 

1959 

1965 

-123.917 

38.587 

9 

23 

4.0 

1933 

1940 

-123.811 

38.478 

10 

23 

1.8 

1623 

1630 

-123.940 

38.393 

11 

23 

24.0 

2191 

2280 

-124.026 

38.493 

12 

23 

21.2 

3111 

3125 

-124.100 

38.591 

13 

23 

18.8 

2443 

2500 

-124.212 

38.700 

14 

23 

16.1 

3197 

3200 

-124.313 

38.808 

15 

23 

12.8 

3257 

3260 


Table 1. Characteristics of the CTD Stations 
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LONG(W) 

LAT(N) 

StaNo 

DATE 

TIME (h) 



-124.456 

38.703 

16 

24 

11.4 

2041 

3500 

-124.332 

38.595 

17 

24 

9.1 

2035 

3445 

-124.239 

38.507 

18 

24 

6.4 

3411 


-124.131 

38.414 

19 

24 


2033 

3430 

-124.027 

38.318 

20 

24 

2.0 

2031 

3430 

-124.149 

38.236 

21 

25 

11.6 

2031 

3600 

-124.248 

38.338 

22 

25 

9.4 

2035 

3670 

-124.356 

38.422 

23 

25 

0.6 

3563 

3650 

-124.467 

38.525 

24 

24 

22.0 

2029 

3605 

-124.574 

38.621 

25 

24 . 

19.6 

2037 

3660 

-124.787 

38.807 

26 

24 

14.5 

2035 

3610 

-124.681 

38.714 

27 

24 

16.9 

2033 

3650 

-124.498 

38.356 

28 

25 

3.8 

2033 

3760 

-124.629 

38.282 

29 

25 

6.0 

2033 

3840 

-124.049 

38.150 

30 

25 

13.6 

2033 

3610 

-123.948 

38.061 

31 

25 

15.7 

2031 

3510 

-123.842 

37.966 

32 

25 

17.9 

2049 

3480 


Table 1. Characteristics of the CTD Stations (Continued) 


The instrument’s pressme, temperature and conductivity sensors were calibrated 
prior to and after the cruise. A General Oceanics rosette sampler was used to collect water 
samples in situ for conductivity calibration after the cruise. The CTD sampling rate was 32 
Hz. The raw data were collected and processed using a software package developed at 
EG&G Marine Instruments. The processing done at this level consisted of removing bad 
data points (outliers) using a first difference criteria to detect suspicious data points and 
linear interpolation to perform the correction. The data were then averaged to every 2 dbar 
and recorded in files for later data analysis. 

An underway data acquisition loop was active during the cruise, automatically 
recording averages of temperature and salinity (among other data) at a depth of 2 m every 
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Figure 2. Position and number of CTD stations. Current meter moorings are 
indicated by the symbol ©. Also shown are SAIL stations (symbol *), the 1000 m 
isobath and the area of study (model domain). 


30 seconds (Ship’s Data Acquisition System, SDAS). These surface data were used to help 
fill data void areas in the interpolation phase, and will be referred to hereafter as SAIL data 
(fi'om Serial Ascii Interface Loop). SAIL data were selected at four positions (SAIL 
stations) separated by approximately the same spacing of the CTD stations, which is about 
14 km (see Figure 2). 

The SDAS also recorded wind speed and wind direction. Figure 3 shows the wind 
conditions during the cruise. Although there was considerable wind variability during the 
survey period, there is no indication that it was in a scale, in space or time, comparable to 
the observed ocean eddy. Thus, although they may be important in the early eddy formation 
stages, the local effects of wind and surface heat flux forcing are neglected in this study. 
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Figure 3. Top: wind vector (direction in degrees true); Bottom: wind speed in m/s 
during May, 1993 off Pt. Arena. The vertical dotted lines mark the beginning and 
end of the cruise. 


A hull-mounted ADCP looking downward recorded velocity data continuously 
following the track of the ship, using accurate navigational information from the Global 
Positioning System (GPS) to compute the absolute velocity. The ADCP is able to measure 
the three components of the current using an array of four acoustic beams; the Doppler shift 
induced in the returning echoes by the moving column of water is processed to compute the 
absolute three dimensional current vector. An appropriate time gating of the echo yields a 
profile of the velocity field in the vertical at several depth bins (for more details on the 
ADCP instrument, see Acoustic Doppler Current Profilers, 1989). 

In addition to the ADCP current observations, moored current meters (conventional 
rotor meters) were also deployed at station numbers 3, 7, 8, 9, 13,17, 18, 19 and 23 (see 
Figure 2). These instruments were deployed at 100,150,300 and 600 m. One mooring also 
had an instrument at 1800 m. Only current meter data from stations 7,8,9 and 13 were 
available for this study. The remaining instruments were recovered at a later time and their 
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data are still being processed elsewhere (Scripps Institution of Oceanography and Oregon 
State University). Table 2 lists the vertical distribution of instruments by mooring, 
identified by CTD station number. 


STATION 

NUMBER 

DEPTH (m) 

100 

150 

300 

600 

1800 

3 


X 

X 



7 

X 

X 

X 

X 


8 

X 

X 

X 

X 

X 

9 

X 

X 

X 



13 

X 

X 

X 

X 



Table 2. Spatial Distribution of Current Meters 


It is assumed that the uncertainty associated with ADCP velocity measurements, 
which can be as large as 10 to 15% of the mean horizontal current (e.g.. Smith and 
Morrison, 1989; Geyer and Signell, 1990), is well beyond the typical magnitude of low 
frequency vertical motions, thus large enough to disregard the use of ADCP vertical 
velocity estimates in this study. However, the observed horizontal velocities from both the 
ADCP and the current meters were used to define the level of known motion (see section 
in.A.2) and to compare model results with observations (see section rV.A.l). Time series 
of temperature and velocity data collected by the current meters were also used in 
estimating the effect of internal tides on the results (section rV.A.3.a). 

B. DATA ANALYSIS 

An inspection of the maximum depth attained on each CTD station in Table 1, 
“Characteristics of the CTD Stations”, shows that whereas the deeper stations can be 
referenced to a common depth of 2000 m, the shallower stations (1 to 5, 8 to 10) do not 
reach that depth due to the shoaling bottom. For reasons that will become clear further on 
(section HI.A), a similarity method was used here to extend those stations’ profiles down 
to 2000 m, to achieve a regular three-dimensional, flat bottom domain. This method 
searches for temperatures in the nearest complete profile equal to the deepest temperature 
in the profile to be extended. After the match is found, the lower part of the former is copied 
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to the latter. The similarity method is attributed to Helland-Hansen, but no citation was 
found. 

The now complete (common depth) CTD data were then analyzed using an optimal 
interpolation method first introduced by Hardy (1971), the multiquadric-biharmonic 
(MQ-B) interpolation method, that has seen several recent applications in meteorology 
(Nuss and Titley, 1994) and oceanography (Chumbinho, 1993). Other authors have applied 
this method to other fields of research, including the one for which the method was initially 
designed, the representation of topography. This method has been studied by several 
authors (Hardy, 1990; Sirayanone, 1988; Nuss, 1994), who believe that it presents 
advantages over other common interpolation schemes (e.g., Barnes and Cressman objective 
analyses) in applications where irregularly spaced, scattered data must be mapped into a 
regular grid. The method behaves very benignly in data void areas, that is, it does not 
introduce unrealistic scales in these areas, as well as it does not introduce scales smaller 
than the spacing between data points. Moreover, the simplest form of MQ-B interpolation 
(the collocation mode, after Hardy, 1990) is of very simple implementation. 

The purpose of using this method in the data analysis was to combine the CTD data 
below the surface with the satellite data at the surface, and to extend the observed gradients 
horizontally outside the CTD area to fill the model domain. However, preference is given 
to the horizontal density gradients because they determine the vertical shear of the currents; 
this limits somewhat the propagation of surface temperatures to deeper levels, which is also 
desirable. The satellite-derived SST are representative of the surface skin depth only, and 
thus may not be related to the temperatures at deeper levels. Haney et al. (1995) show, in a 
study that attempted to recover temperature profiles from shallow CTD profiles using 
empirical orthogonal functions, that CTD casts to shallow depths less than 100 m are not 
capable of representing the whole water column. Thus, the main result of this MQ-B 
analysis is the correct representation of the horizontal density gradients determining the 
vertical shear. 

1. Multiquadric-biharmonic interpolation method 

This section gives a description of the MQ-B method, following closely that 
provided by Nuss (1993). Consider a three-dimensional irregular surface containing a field 
sampled at N observation points. For each observation point we know the location, (x,y,z), 
and the value of the field, H. This irregular surface can represent any scalar (and, in more 
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recent applications, also vector) field. The purpose is to reconstruct the irregular surface 
based on the observed data. 

The interpolation of the field value in any arbitrary point X is based on a weighted 
sum of radial basis functions: 

N 

H(X) = (2.1) 

/ = 1 

where a, is the weight of the kernel function Q whose argument is the radial distance 
between the field point X and the observation point X,-. In the multiquadric-biharmonic 
interpolation method, the kernel function is a quadric function: 

Q{X-X^) = [(X-xf + c^]^^. (2.2) 

The constant c is a small arbitrary number to keep the function from vanishing at 
the observation points. For the three-dimensional interpolation case, the kernel function is 
therefore: 

1/2 

QiX-Xi) = [(x-xf+iy-yf+iz-zf + c'^] . (2.3) 

To solve for the weights, we apply Equation 2.1 to all N observation points, 
resulting in a set of N equations with N unknowns: 

N 

HiXj) = Y^iQi(Xj,yj,Zj). (2.4) 

1=1 

In matrix notation, the solution is given by: 

a = Q"^ • Hp (2.5) 

where Q is the square matrix of equation coefficients and Hf is a vector of the observed 
field values. This matrix can be inverted using any of the mathematical routines available 
for solutions of linear systems of equations (e.g., Linpack, 1978, or Matlab, 1994), 

Once the weights are known, the interpolated solution to any point in the field (in 
usual applications, regularly spaced grid points) is then given in matrix notation by: 

Hg = Qg • a. (2.6) 

where: 
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(2.7) 


1/2 

Qg = [(Xg-xif+ (yg-yif+ (z^-z,)^ + c^] , 

(Xg,yg,Zg) being the coordinates of each desired field point. The matrix on the left hand side 
of Equation 2.7 has as many rows as there are field (grid) points, but the number of 
columns is always equal to the number of observation points N. In other words, the 
solution depends only on the number of observation points available and is independent of 
the size of the grid to be interpolated. The two-dimensional problem is just a particular 
case of the more general three-dimensional case, with all z coordinates set to zero. 

The condition that the sum of the weights be equal to zero (Hardy, 1990) is 
implemented by adding an extra row and column to Q as follows: 



[ 1 - 1 ] « 


( 2 . 8 ) 


where the vectors have length N, and by appending a zero to the data vector Hf. This, in 
turn, introduces a new term in the weight vector: 



^2 





(2.9) 


The new term represents the average value of the entire observation data and thus 
must be added to the final interpolated field. The multiquadric matrix in Equation 2.6 
maintains the dimensions mentioned above (i.e., MxN, where M is the total number of grid 
points), since for the field interpolation only the first N weights are used. The MQ-B 
interpolation equations are: 

Hg = Qg«a + ao, (2.10) 

with i - 1,..., N and g - 1,..., M. 

The method presented above is valid for noise-free observations. Unfortunately, 
real data include observational errors that must be taken into account. This is achieved 
through the introduction of an uncertainty parameter and a smoothing parameter (Nuss and 
Titley, 1993) in Equation 2.4. In matrix notation: 

Hj = [Q + • 5.^.] • a , (2.11) 
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where X is the smoothing parameter, is a vector of uncertainties in the observations 
(mean-squared observation error) and 6y is the Kronecker delta function. The difference 
now is that the diagonal elements of Q are modified to include a factor involving those 
two parameters. Note that the uncertainties vector can have different values for different 
observation points, a useful property that gives the user of the method enough freedom to 
optimally set the parameters for best fitting any part of the grid. 

2. Data interpolation 

The multiquadric-biharmonic interpolation method described above was applied to 
the temperature and salinity data to obtain an analyzed three-dimensional density field in 
the model domain. Chumbinho (1993) presents some thoughts about the MQ-B method, 
namely the choice of the interpolation parameters c and X, and the choice of an appropriate 
vertical scale exaggeration. The interpolation parameters were set to values close to those 
used in other applications of MQ-B (Nuss, 1994; Nuss and Titley, 1994), in a way that 
guarantees the best representation of local maxima and minima, without distortion of the 
interpolated field in regions where data is scarce or not available. The vertical exaggeration 
is necessary to emphasize the correct representation of the horizontal gradients over the 
vertical ones, which are usually much larger. 

The number of temperature observations is much larger than the number of salinity 
observations because of the use of satellite-derived sea surface temperatures. As noted 
previously, these were used in an attempt to compensate for the lack of CTD data in several 
parts of the domain. For this reason, the data interpolation was carried out differently for 
these two state variables. 

a. Temperature interpolation 

The CTD temperature data were first decimated along each profile to the 
same 27 depth levels used by the circulation model (-3.5, -12.5, -24.3, -38.0, -54.9, -75.9, 
-101, -131, -166, -206, -252, -303, -360, -423, -493, -570, -653, -744, -842, -948, -1060, - 
1182, -1312, -1450, -1596, -1750 and -1916 m; also, see section ni.A). This decimation 
was performed following the filtering criteria defined by Crochiere and Rabiner (1983) for 
non-uniform decimation, and according to the idea that the MQ-B interpolation works best 
when the data points coincide with grid points (Hardy, 1990, p. 164). This resulted in 32 
in situ temperature data points at each of the 27 vertical levels, for a total of 864 points. 
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The satellite-derived sea surface temperature data (SST) were sampled 
directly from the multispectral satellite image equivalent to the single-channel image in 
Figure 1, at a resolution of 48x39 points, respectively in the alongshore and cross-shore 
directions, covering the model domain represented by the rectangle in Figure 2, using a 
geographical to Universal Transverse Mercator (UTM) coordinates transformation and 
posterior rotation. The origin of the UTM coordinates is the northwestern comer (or lower 
left comer, in the figures). The total number of SST data points is thus 1872, which is much 
larger than the 32 temperature data points at the second level. 

The MQ-B performs best when the data are irregularly spaced, without 
sharp contrasts in spatial data density (Hardy, 1990). This was readily observed in the first 
few interpolation mns without any pre-processing: the three-dimensionally interpolated 
temperature field presented severe distortions in the upper levels, caused by the contrast in 
the number of total data points at the surface and at the level immediately below. This 
problem lead to the necessity of pre-processing the temperature observations. 

The SST data were thus two-dimensionally decimated by a factor of seven, 
resulting in a rectangular array of 8x7 data points covering the same area, again following 
the criteria defined by Crochiere and Rabiner (1983); at the same time, a two-dimensional 
MQ-B interpolation of temperature was carried out at each level, using very small 
parameters (c=lxl0''^, X infinitesimally small), to artificially increase the number of 
observations at depth from 25 to 64 (see below). The interpolation parameters mentioned 
above introduced no smoothing in the interpolated field, and caused the resulting field to 
follow the data points exactly (the interpolation results were checked using three- 
dimensional visualization of the data points against the interpolated field). Only the 
decimated temperature profiles corresponding to stations 1 to 25 were used in this 
two-dimensional interpolation, and the resulting field was mapped into a square 8x8 grid 
about the same size as the horizontal area occupied by these stations, as illustrated in 
Figure 4. 

The last step in the temperature interpolation was to combine the 8x7 
decimated SST, the 8x8x27 two-dimensionally interpolated temperatures, and the 7x27 
points from the remaining seven CTD stations (numbers 26 to 32) in a three-dimensional 
MQ-B interpolation, using c-2xl0 ^ X=lxl0'^ and an aspect ratio of about 1/60 (achieved 
by scaling the horizontal dimensions to unity and leaving the vertical dimension unaltered). 
The resulting three-dimensionally analyzed temperature field followed the data points 
closely; profiles of analyzed data were compared with nearby CTD profiles, showing a 
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Distance Alongshore (km) 


Figure 4. Schematic representation of the 2-D temperature interpolation at each 
level. The dots are the model domain grid points, the crosses are the CTD stations 
and the rectangle is the area covered in the 2-D interpolation, mapped by a 8x8 
grid. 

good match. The maximum differences between analyzed and observed global maximum 
and minimum temperatures were 0.11 and 0.07 degrees Celsius, respectively. This 
analyzed field had 48x39 points at each of the 27 levels. 

b. Salinity interpolation 

The MQ-B interpolation of salinity was considerably simpler than the 
interpolation of temperature, because there were no sharp contrasts in the number of 
observations between levels. The CTD salinity profiles were non-uniformly decimated in 
the vertical, in a manner similar to the temperature vertical decimation, to the same 27 
depth levels, resulting in 32x27 observation points. Four other salinity observation points 
were included in the salinity data vector, corresponding to the SAIL stations, and surface 
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salinities were taken from an average of the first two CTD records on each profile. The total 
number of salinity observations was, therefore, 900 data points. 

These observations were three-dimensionally MQ-B interpolated using 
c=lxl0'^ ^=lxl0‘^ and an aspect ratio of about 1/60, similar to the temperature 
interpolation. The slightly different value for c yielded better results, maybe because of the 
different spatial distribution of salinity observation points compared to that of the 
temperature observations. The resulting three-dimensionally analyzed salinity field also 
followed the data points closely, with very good match between analyzed profiles and 
nearby CTD salinity profiles. The maximum differences between analyzed and observed 
global maximum and minimum salinities were 0.01 and 0.04 psu, respectively. The 
dimensions of this field were the same of the temperature field (48x39x27). 

3. Analyzed three-dimensional density field 

Once the temperature and salinity data, from the various sources described above, 
were interpolated and analyzed, the three-dimensional (48x39x27) density field was 
computed using the formulae for p(S,T,p) in the International Equation of State (lES 80) 
found in Gill (1982). The values of density were then converted to sigma values: 

o(s,r,p) = [p(s,r,/?)-iooo] (2.12) 

where S is salinity, T is temperature and p is pressure. 

The result was a three-dimensional, analyzed density field from which the three- 
dimensional currents were diagnosed using the primitive equation model, and the 
quasigeostrophic vertical motion was computed from the quasigeostrophic omega 
equation. An implicit assumption of synopticity in the data was made throughout the 
analysis. The hydrographic data were collected during a period of three days, during which 
several processes can act to distort the sampled field. Among these are the rotation and 
translation of the nearly circular eddy (Figure 1), and the propagation of internal waves in 
the domain sampled by the CTD. The most important internal waves that could be present 
are internal tides, which have been known to contaminate data sets in other areas. The 
possible effects of these processes are analyzed and discussed in section rV.A.3.a. 
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III. METHODS 


After the data were analyzed and the three-dimensional density field computed, a 
field of dynamically balanced currents in the area of study was obtained by initializing a 
primitive equation model. To help interpret the results, the vertical velocity was also 
diagnosed from the (quasigeostrophic) omega equation. Both the primitive equation model 
and the quasigeostrophic omega equation were applied in the volume defined by the 
analyzed density field (48x39x27, respectively alongshore, cross-shore and in the vertical) 
whose horizontal extension is depicted in Figure 1 (hereafter domain). This domain was 
chosen to cover a larger horizontal area than the domain defined by the CTD stations to 
allow for model boundary effects and possible eddy translation. 

This chapter describes several aspects relevant to each of these solution methods. 

A. CIRCULATION MODEL 

The objectively analyzed density field was used in a circulation model for the 
diagnostic computation and analysis of the dynamically balanced, three-dimensional 
velocity and density fields corresponding to the slow mode assumed captured in the 
observations. 

1. Description 

The model used in this study is the multi-level, primitive equation model presented 
by Haney (1974); it utilizes a B grid in the horizontal and sigma coordinates in the vertical, 
with vertical staggering of the vertical component of velocity. The inviscid, hydrostatic and 
Boussinesq assumptions are used to simplify the set of governing equations. This means 
that the ocean is considered incompressible, and density is replaced by a constant except 
where it is multiplied by the acceleration of gravity. 

This model has been improved with the modifications introduced in the convective 
adjustment algorithm by Adamec et al. (1981) and in the numerical computation of the 
pressure gradient force in sigma coordinates and in the vertical discretization scheme by 
Haney (1991). The convective adjustment scheme imposes dynamic stability by using an 
algorithm that conserves heat and momentum and mixes these two properties between 
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levels at the same ratio, requiring a final Richardson number (after the scheme has been 
applied) greater than the critical value of 0.25. The new vertical discretization and pressure 
gradient force computation were designed to improve this computation over steep 
topography. 

In this study, the dynamics of the observed eddy will be studied under the working 
hypothesis that surface forcing (wind stress and buoyancy fluxes) and variations in bottom 
topography can be neglected. This simplifies the model by permitting the use of a single 
thermodynamic variable, o, given in Equation 2.12, and a flat bottom. The possible effects 
of these assumptions will be discussed in Chapter V. 

2. Initial conditions 

To achieve dynamic initialization of the primitive equation model, information on 
the initial mass and velocity fields is needed. The former is the objectively analyzed 
density, or a; the latter can be set in one of two ways, if no direct current observations are 
available in the model domain. The initial shear velocity field can either be set to zero, 
implying that no a priori knowledge of the initial velocities is available, or it can be 
specified from the initial mass field using the hydrostatic and geostrophic relations. These 
two initial velocity fields have different implications in the final balanced slow mode 
solution, as will be seen in section IV.A.2. 

The depth averaged flow is also computed from the hydrostatic and geostrophic 
relations, assuming a level of known motion determined by the current meter data. The 
resulting vertically averaged flow, being nondivergent, is then used to compute a 
streamfunction for the vertically averaged velocities. This streamfunction is then kept 
constant during the short model initialization run. 

The level of known motion was chosen to be 2000 m, with a current of 2 cm/s 
directed alongshore to the northwest, based on current meter data. This apparently simple 
method was preferred over a more complex method of determining a level of known motion 
by using the horizontal ADCP currents. This approach would have required an objective 
analysis of the ADCP currents and the solution of a Poisson equation for the rotational part 
of the velocity field with consequent specification of boundary conditions for the 
streamfunction. The correct determination of the depth averaged flow is important because 
of its effect in the advection of relative vorticity perturbations; however, it will be shown 
that the solution is barely sensitive to variations in the setting of the bottom current. 
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3. Boundary conditions 


The model has open lateral boundaries governed by modified Ross and Orlanski 
(1982) radiation conditions, hereafter RO boundary conditions. This modified scheme 
serves to radiate away from the domain any outgoing waves (mainly gravity waves), and to 
include inward advection of properties specified at the boundaries. The RO scheme can be 
summarized in three steps: 

• computation of a normalized phase velocity 6 for propagating waves (see 
Equation 3.5 below); 

• for inflow points (5 < 0 ), the vorticity, divergence and density are set to the 
observed boundary values; 

• for outflow points (5 > 0 ), the vorticity, divergence and density are advected from 
the interior using equations of the type: 


dt ^dn 



(3.1) 


In Equation 3.1, | is either divergence, vorticity or density, and the spatial 
derivative is taken along an outward normal n to the boundary. Once the values of 
divergence and vorticity are known at the boundaries, their definitions are used to compute 
the normal and tangential components of the baroclinic horizontal velocity. For outflow, 
density is also advected along the boundary with an equation similar to Equation 3.1, with 
the spatial derivative generalized to include its tangential counterpart, and a new phase 
velocity 82 computed as the spatial four-way average of the surrounding tangential 
velocities at each boundary point. 

The boundary equations for vorticity and divergence in the outflow case were 
modified from the original equations given in RO, since their equations were mainly 
designed to deal with gravity waves, with divergence waves propagating faster than 
vorticity waves. RO use the phase speed 5 to advect vorticity, and Ax/At to advect 
divergence, where Ax is the grid spacing and 1st is the time increment. This means that the 
divergence at the boundary is set equal to the value at the adjacent interior point, computed 
from the prognostic equations governing the interior. 

The velocity field in the model initialization runs described below presents 
noticeable inertial oscillations, where divergence and vorticity waves propagate at 
comparable speeds. The different advection speeds for one and the other in the RO method 
caused the components of horizontal velocity to become excessively large on the 
boundaries. Because of this it was necessary to include the inertial terms and a damping 
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The computation of an appropriate phase speed 8 was also revised from RO. 
Following the work of Durran et al. (1993), the computation of the phase speed was 
changed to: 

8 = CJ, + V • n (3.5) 


where V is the velocity vector at the boundary, n is an outward unit vector, and Q is a 
constant, set as to guarantee propagation of waves out of the domain during most of the 
integration time (equivalent to the shallow water phase speed of gravity waves). V • n is a 
four-way average of the normal velocity; for instance, for the western boundary V • n is 
computed as: 


V • n = -i -h u^^]_ 1 -h u[^] -h 1 ] 


The superscript (n) means values of the variables at time n. This makes the 
advection-type Equation 3.1 a completely forward-in-time scheme (as it is), which is 
acceptable for upstream advections (it is used to compute unknown points at the boundary 
from known interior points, in outflow situations). The final vorticity and divergence 
equations are, respectively (in finite differences form): 


- 8 !|| - 4 ">] 

fs"” \.D^U +/A(Ci"+‘>-7|A(|DW 
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where the superscripts (n) and (n+1) mean time level n and time level n+i, respectively, 
and the subscripts B and B+1 mean the boundary point and the adjacent interior point. 

Note that because Equation 3,7 uses implicit time differencing for the Coriolis 
terms, it has to be solved simultaneously for the new time levels. This is done at each time 
step along all four boundaries, except in the comers. The vorticity and divergence in the 
comers are set equal to the average of the two adjacent points that have already been 
computed by advection normal to their boundaries. 

For use in the case of inflow at the lateral boundaries, the vorticity was specified 
from its initial value and the divergence (also from its initial value) was zero. These 
conditions were held constant during the short model initialization mns. 

The boundary conditions for the vertical velocity in this study were zero at the 
bottom and at the surface. 

4. Digital filter initialization 

The objectively analyzed density field and geostrophically balanced currents were 
used as initial conditions in the model. Separate initialization runs were also made with the 
initial currents at rest. These fields are not in balance, or at least not in the state of near 
balance that characterizes a quasigeostrophic ocean, A digital filter initialization (DFI) 
procedure introduced by Lynch and Huang (1992) was used to arrive at a dynamically 
balanced field of density and currents. The DFI procedure, which removes high-frequency 
oscillations from the model fields, consists of performing two short (12 h), adiabatic 
(physics turned off), forward and backward integrations of the equations, starting with the 
same initial conditions. The resulting model variables at each time step from t =-12 hio 
t =12 h are multiplied by a spatially uniform, time varying weight. The time series of 
weights corresponds to the coefficients of a digital lowpass filter effectively being applied 
to the computed fields, i.e.: 

N 

qin) = '^q{k)h(n-k) (3.8) 

k = -N 

applied only for n equal to zero. In the equation above, q is any of the variables m,v or o, 
while h(n) is the sequence of filter coefficients, and N is the number of time steps (equal to 
the integration time forward or backward divided by the time step, which was 300 
seconds). This procedure is illustrated in Figure 5. The figure shows an example of a time 
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Figure 5. Example of the DFI procedure applied to a typical vertical velocity time 
series. The solid line is vertical velocity in m/day; the dash-dot line is the sequence 
of filter coefficients (scale on the right) multiplying that time series. 


series of q, in this case the vertical velocity, at one of the grid points, and the time filter 
that is to be applied to it. The output from Equation 3.8 is q (0), the filtered vertical 
velocity at time zero at that grid point. 

The digital filter utilized by Lynch and Huang (1992) consists of a simple lowpass 
prototype multiplied by a Lanczos window to reduce truncation effects. The same authors 
have already proposed other filter design techniques (Huang and Lynch, 1993). The 
possibilities in this area are broad, since various filter responses and designs can be tested, 
as demonstrated by Viudez et al. (1994). In the present study, the filter was implemented 
as in Lynch and Huang (1992): 


h{ri) = 


sin[n 7 i/(iV+l)] sin(ne^) 
nn/{N+\) nn 


(3.9) 
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where 0^ is the cutoff digital frequency. When this filter is used with a cutoff period equal 
to the total time spanned by the two integrations (24 hours) and a time step of 300 seconds, 
the magnitude of the frequency response shown in Figure 6 is obtained. It is seen that the 



Figure 6. Frequency response of a lowpass filter with cutoff period 24 h, sampling 
period 300 s (Lanczos window). The dotted line marks the digital cutoff 
frequency. The dash-dot line marks the response for a period of 12 h. 


cutoff period is damped by about 5 dB. It is possible to design more efficient filters that 
would achieve better frequency responses with fewer coefficients. As explained above in 
Equation 3.8, the number of coefficients in the filter is directly related to the length of the 
forward and backward integrations. Thus, fewer coefficients means less computer time in 
the initialization process. On the other hand, it also means less time for ageostrophic spin- 
up (development of secondary circulations from the unbalanced mass and velocity fields; 
see section IV.A.2). A more efficient filter design (for example, using Butterworth or 
elliptic lowpass filters) would have to take this constraint into account in the computation 
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of the coefficients to allow enough integration time for spin-up (usually no less than 3 to 6 
hours, setting a minimum limit of 36 to 72 coefficients for a finite impulse response filter). 

The application of the filter to the model output does not need to be done within the 
circulation model, as long as enough computer storage space is available to save the model 
output at each time step. In addition to using different filters, the idea of using the 
dynamically balanced (filtered) mass and velocity fields as initial conditions for a second 
DFT (integration-filter) procedure is being tested and verified with model generated data, in 
a separate study. 

In the original paper. Lynch and Huang (1992) show that the DPI procedure is 
equivalent to a nonlinear normal mode initialization, with the advantage of having a much 
simpler implementation. The normal mode initialization compares the magnitude of the 
modes resolvable by a model through long-term integrations (gravity wave modes versus 
Rossby modes); the same modes are fitted to the analyzed data, and the components that 
have negligible contributions to the long-term solution (i.e., the gravity waves) are 
removed from the data. The nonlinear variant of the normal mode initialization takes care 
of the re-generation of high frequencies that takes place even after those frequencies have 
been removed from the data, because of the nonlinear terms in the model equations. 
Therefore, the model is initialized with the slowly varying, almost nondivergent basic state 
captured in the data. The DFT procedure removes the high frequency modes with the 
application of the filter. The lowpass filtered fields then represent the dynamically balanced 
slow mode at the time of the observations. 

Horizontal and vertical diffusion was turned off to enable the short backward and 
forward integrations centered at time zero. When compared to advective processes in 
regions of strong mean flows, such diffusion is negligible, as noted by Schmitz and Vastano 
(1975) in a numerical study of Gulf Stream rings over time scales of two months. 

B. QUASIGEOSTROPHIC OMEGA EQUATION 

The second method of dynamical analysis, which was used to diagnose only the 
vertical velocity from the data, was to solve the quasigeostrophic (QG) omega equation. 
The assumptions allowing the simplification of the complete set of equations governing 
geophysical fluid dynamics into QG equations on a /-plane are the following: 

• Hydrostatic balance; 

• Geostrophic velocity and vorticity are used everywhere except in the divergence 
term in the vorticity equation; 
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• Adiabatic flow; 

• Neglect of the vertical advection of vorticity and the tilting of vorticity; 

• Relative vorticity much less than planetary vorticity in the divergence term; 

• Planetary vorticity is constant; 

• Static stability is a function of z (depth) only. 

An important concept in QG theory is that vertical motion is entirely due to the 
ageostrophic part of the currents, because the geostrophic currents are nondivergent. Thus, 
using the assumptions above, information on the vertical motion is contained in the 
vorticity equation relating local changes to vorticity advection by the geostrophic flow; 


where is the vertical component of the geostrophic relative vorticity defined as: 




dx 




(3.11) 


(the last equality follows from the curl of the geostrophic balance), g being the 
acceleration of gravity and/the planetary vorticity (2Qsin(t)). 

Equation 3.10 has been used to estimate vertical motion after substituting the 
continuity equation for incompressible fluids 

v.y = -^, (3.12) 

oz 

for the horizontal divergence. The resulting equation for w requires knowledge of the local 
variation of vorticity with time, which is not always available. A more useful equation 
without any time derivatives (diagnostic equation) can be derived by combining Equation 
3.10 with another equation relating vertical velocity to hydrostatic pressure perturbations, 
the buoyancy equation (Gill, 1982, page 529): 

Po^,v + [i+v,.v]|=o (3.13) 


where is the buoyancy frequency defined as: 


j2 _ 8 


Pb^z’ 


where p is pressure, p is density, Po is a reference density, w is the vertical velocity, and \g 
is the geostrophic velocity vector. The Laplacian and del operators are two-dimensional 



(horizontal). Equation 3.13 is the nonlinear QG version of the hydrostatic buoyancy 
equation. 

With the help of the hydrostatic relation 

I = -Pg (3.15) 

the buoyancy equation can be written as: 

• V p (3.16) 

From the thermal wind relations 


dz P(/dy 

^ _x?£ 

dz p(/dx 


(3.17) 


an expression for the vertical derivative of geostrophic relative vorticity can be derived: 


- 

dz ~ 


— §-V^p 

Pof ^ 


(3.18) 


Taking the z-derivative of Equation 3.10, the horizontal Laplacian of Equation 
3.16, using the results in Equations 3.12 and 3.18, rearranging terms to eliminate the time 
derivative and using the QG assumptions of constant planetary vorticity (fo) and 
horizontally constant stability, the quasigeostrophic equation for vertical velocity is 
obtained (omega-equation): 

'4^+A('v3jH.=/„^(v^.vy + ^vi“(v^.vp) . (3.19) 


or, in terms of the geostrophic vorticity defined in Equation 3.11: 



Equation 3.20 shows that vertical motion can be diagnosed (no time derivatives) 
from the instantaneous pressure field only, since this variable determines both the 
geostrophic currents and the geostrophic relative vorticity in the CJG framework. The term 
on the left hand side contains only second derivatives of vertical velocity; for disturbances 
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with wave-like character, it is proportional to -w. The contribution of each of the terms on 
the right hand side of this equation to the vertical motion will be discussed in section IV.B. 

The data available for this study is a snapshot of the density field, with considerable 
detail in both the vertical and the horizontal (after the analysis described in section II.B). 
This is an ideal situation for the application of Equation 3.19 for the computation of vertical 
motion. However, it presents some problems because the forcing terms on the right hand 
side (differential vorticity advection and Laplacian of thickness advection) sometimes tend 
to cancel each other, leading to errors in the numerical procedure to compute the response 
on the left hand side. A more stable alternative, the Q-vector form, was introduced by 
Hoskins et al. (1978): 

(iV^V2+y2_^)H- = 2V«Q, (3.21) 

dz 

where 


Q 


Po 




(3.22) 


This form of the (o-equation (Equation 3.21 above is the so-called Trenberth form, 
where deformation terms arising from taking the Laplacian of a product have been 
neglected) shows that a convergent Q-vector field implies upward vertical motion. 
Moreover, the vectors themselves are parallel to the ageostrophic velocities responsible for 
the vertical motion. This Q-vector form presents additional advantages since the forcing 
term is exact within the (^G context. Neglecting the deformation terms could be a source of 
error in the vicinity of strong jets or frontal systems (Carlson, 1991), however that is not a 
problem in the present case because even a coastal jet has sufficiently small velocity 
magnitudes in an eastern boundary current. The forcing functions can also be evaluated on 
a single level surface. 

Several aspects of the QG assumptions can cause the results from the (o-equation 
to be either under- or over-estimated. The departures from QG arise in regions where: 
1) diabatic effects are strong; 2) static stability is not constant in the horizontal; 
and 3) ageostrophic motions due to curvature are so strong that neglecting ageostrophic 
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advection becomes an important source of error. The relative importance of these terms will 
be considered in the discussion of the QG results (section IV.B). 

The o)-equation was solved using both the original MQ-B analyzed density field, 
and the dynamically balanced density field resulting from the DPI procedure of section 
in.A.4, to compute the forcing (right hand side of Equation 3.21). In both solutions the 
condition of w=0 on all boundaries was used. 
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Figure 7. Horizontal cross section of the dynamically balanced horizontal velocity 
field at a depth of 100 m, in cm/s (red arrows). The longest arrow is about 50 cm/s. 
The symbol © represents the CTD stations. The yellow arrows are the cuixent vectors 
from the current meters at the same depth (stations 7, 8, 9 and 13). The white line is a 
streamline of the nondivergent component of the velocity. The dashed lines are 
included for reference and represent the vertical cross-sections shown in other figures. 
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The mean flow is not circular around an hypothetical eddy center; instead, it has meanders, 
with at least five well defined troughs. These smaller scale features (hereafter sub-eddy 
scales) have been found in other systems, particularly in frontal systems (e.g., Fiekas et al, 
1994) where submesoscale circulations are part of the mesoscale meandering of the front 
itself. This finding is partly corroborated by an inspection of Figure 1. This figure shows 
several feeder bands of water being entrained into the eddy circulation, suggesting that the 
pattern of motion around the eddy is not completely circular. 

A comparison was made between the DFI solution of the horizontal velocities and 
the currents observed by the hull mounted ADCP (section n.A.2). Figure 8 shows a 
horizontal map of the dynamically balanced horizontal currents and the horizontal ADCP 
vectors for the bin centered at 100 m depth. In general, the DFI solution and the ADCP 
currents compare well in magnitude and direction. This is especially so in the northeast, 
northwest and southwest parts of the domain, even though the DFI solution in the 
northeastern part of the domain (upper side in Figure 8) presents some undesirable effects 
from the flat bottom assumption (see below). The region to the southeast of the eddy is an 
area where practically no CTD observations were made, and thus the solution is almost 
entirely dominated by the performance of the MQ-B method in data void areas, which was 
said earlier to be very benign. This means that the density gradients observed by the CTD 
stations at the southeast edge of the eddy become smaller and smaller as a function of 
distance from the eddy, and the general circulation will therefore tend to remain weakly 
cyclonic. An exception is the small anticyclonic feature centered at (31,15), which is 
persistent in the DFI results down to several hundred meters and was never sampled by the 
ADCP (results not shown). Since this position coincides with a CTD station, the resulting 
anticyclonic circulation is the response to the observed density gradients between stations 
22, 21 and 30, extending quite deep. 

In the analysis of Figure 8 it should also be kept in mind that the ADCP observes 
the total velocity with a measurement uncertainty that is highly variable and, in most cases, 
of the same order of magnitude of the signal that is being computed in the DFI solution (the 
ageostrophic part of the flow). Also, the DFI solution is representative of the slow mode at 
the time of the observations, whereas the ADCP is observing all frequencies. The ADCP 
data plotted in Figure 8 was spatially averaged, resulting in effective time averages that 
range from 15 minutes to about 2 hours. 

The meanders in the mean flow have a distinct field of vertical motion associated 
with them. This is visible in Figure 9 which shows the vertical velocity at 100 m (colored 
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Figure 8. Florizontal cross section of the dynamically balanced horizontal velocity 
field at a depth of 100 m, in cm/s (red arrows), with superimposed ADCP horizontal 
velocities for the depth bin centered at 100 m (yellow arrows). 
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Figure 9. Horizontal cross section of the dynamically balanced vertical velocity field 
at a depth of 100 m, in m/day (colored field), with superimposed horizontal velocities 
at the same level. The symbol © represents the CTD stations. The red arrows are the 
current vectors from the current meters at the same depth (stations 7, 8, 9 and 13). The 
black bold line is a streamline of the nondivergent component of the velocity. The 
black thin lines are included for reference and represent the vertical cross-sections 
shown in other figures. 







field) with the horizontal velocity (black arrows) and the same streamline as in Figure 7 
superimposed. The w field is distributed mostly around the eddy in patches of radially 
aligned, alternating upwelling and downwelling centers, with maximum vertical velocities 
at this level of about ± 20 m/day. The average vertical velocity within a radius of 35 km 
from the center of the eddy, however, is in the order of 0.1 m/day. 

These sub-eddy scale patches of rising and sinking motion are related to meanders 
in the otherwise circular eddy and can be explained by QG theory by considering the 
changes of vorticity in a water parcel following the fluid through the meanders (Holton, 
1979). Figure 10 shows a typical meandering system with cells of negative and positive 

A C 



Figure 10. Meandering jet with alternating patterns of negative (regions A and C) 
and positive (region B) vorticity. The changes of vorticity going from A to B to C 
cause downwelling upstream of the trough and upwelling downstream of the trough. 

vorticity. In region A, the parcel acquires negative relative vorticity; as it moves to region 
B, it encounters regions of increasing relative vorticity. The QG vorticity equation 
(Equation 3.10) can be written: 

f = -/V.V = /|. (4.2) 

Hence, from A to B the left hand side of Equation 4.2 is positive; with zero boundary 
condition at the surface, the right hand side shows that w is negative at depth, that is, from 
A to B there is sinking motion. If the meander is not growing or decaying, the zero vertical 
velocity point will be exactly at the trough, because in this case the left hand side is zero 
there. 




The opposite happens as the parcel moves from region B to region C; the change in 
the relative vorticity of the parcel is now from positive values to negative values, causing 
the left hand side of Equation 4.2 to be negative. Using the same boundary condition on the 
right hand side a positive value for the vertical velocity is now found at depth. Thus an 
upwelling cell lies between region B and region C. The flow is thus three-dimensional, with 
the particle describing not only the meanders in the mean flow but also some vertical 
excursion in both senses (up and down). In this simple but illustrative model the vertical 
displacement of the particle is: 

® inversely proportional to the magnitude of the horizontal mean flow between A 
and B or between B and C; 

» directly proportional to the magnitude of the vertical velocity. 

The second factor is proportional to the rate of change of vorticity following the 
motion and the depth over which the vorticity change occurs. All else being equal, large 
vertical scales tend to have stronger vertical velocities. In the Gulf Stream meanders, for 
instance, vertical motions of 100 to 200 m/day have been found (Lindstrom and Watts, 
1994). In this region the vertical motion may be enhanced by the large density gradients 
across the meanders, which cause the changes in vorticity to take place over appreciable 
changes in water column depth, even for small displacements from the jet axis. The rate of 
change of vorticity is largest in regions of strong changes in curvature, i.e., for smaller 
horizontal scales. 

The upwelling and downwelling patterns described above can be identified in 
Figure 9. The patches of vertical motion are correlated with the meanders of the mean flow 
depicted by the streamline, with downwelling upstream and upwelling downstream of the 
troughs. Similar patterns of vertical motion have been identified in the meanders of the Gulf 
Stream (Leach, 1987; Bower and Rossby, 1989; Bower, 1989). The relative magnitudes of 
the horizontal and vertical velocities also show that the vertical excursion of the water 
parcels in these patches is small. A typical horizontal length scale for the patches is about 
20-25 km. This distance is covered in about 15 hours at a speed of 40 cm/s; in 15 hours a 
parcel moving at this speed will be displaced a maximum of only 12 m in the vertical. 

Another important aspect of the description of this eddy is its vertical extension. 
This can be appreciated in Figures 11 and 12. The former is a vertical cross-section of the 
vertical velocity along the x-axis (alongshore), corresponding to the horizontal line in 
Figure 7; the latter is a cross-section of the vertical velocity along the y-axis (across shore), 
corresponding to the vertical line in Figure 7. Both have units of m/day, although the color 


37 






Figure 11. Vertical cross section of the dynamically balanced vertical velocity 
(colored field) in the alongshore direction (from the northwest (left) to the southeast 
(right)), in m/day, with depths in meters. Also shown are the two-dimensional (m,w) 
vectors in the same section (black arrows). Note that the vertical axis is not linear. 







Figure 12. Vertical cross section of the dynamically balanced vertical velocity 
(colored field) in the across shore direction (offshore to the left), in m/day, with depths 
in meters. Also shown are the two-dimensional (v,w) vectors in the same section (black 
arrows). Note that the vertical axis is not linear. 
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scale is not the same in the two figures; also, note that the vertical axis is not linear (the top 
500 m take about half the axis). The vertical extent of the eddy is shown to be about 500 m 
(extending in some points to almost 1000 m), with maximum vertical velocities occurring 
between 50 and 150 m depth. The patches of rising and sinking motion are not uniquely 
related to divergence and convergence of the along-section vectors superimposed on the 
figures because there may be compensation in the across section flow. However it can be 
seen that the strong local maxima of rising motion generally have associated with them 
divergence of the velocity vectors above and at the same level while the local minima are 
generally associated with convergence above and at the same level. 

It is also noticeable that in Figure 12 there are large horizontal vectors in the far 
right of the figure, that is, the boundary of the model closer to shore. This northeastern 
boundary, and the northwestern boundary, are regions of strong tangential velocities, at 
least in the upper levels as can be seen in Figures 1 and 7. This causes some problems at 
the boundary, which use an advective velocity based only on the normal flow in the 
boundary condition (section ni.A.3). Furthermore, this is also the region where the 
Helland-Hansen method was used to extend the inshore CTD profiles to a common 
maximum depth (section n.B), thus introducing horizontal density gradients that may not 
be real. These horizontal gradients cause unrealistic vertical shear in the horizontal 
velocities through the thermal wind relations. For these reasons, the model results in the 
northeastern part of the domain and closer to the shore may not be reliable, and the 
discussion is therefore limited to regions where the joint effect of boundary and missing 
topography is clearly nonexistent or minimized. 

2. Geostrophic adjustment 

It was mentioned in the preceding section that the results described above were 
obtained with the DFI procedure in a model whose initial baroclinic velocity field is 
geostrophic. Another possibility is to set the baroclinic currents to zero at the initial time. 
This different initial condition will also yield a final slow mode; the choice between one 
and the other is related to the process of geostrophic adjustment itself. 

Geostrophic adjustment is the process by which an initial disturbance in a rotating 
fluid releases part of its potential energy into kinetic energy, until a steady equilibrium is 
achieved (Gill, 1982). The real ocean is in a constant process of geostrophic adjustment due 
to multiple forcing. In a simulated ocean without forcing, the geostrophic adjustment takes 


40 



place during the spin-up phase, where the imbalance between the initial density and mass 
fields is geostrophically adjusted. During the process, about one third of the released 
potential energy is converted to kinetic energy of the equilibrium solution, and the 
remaining two thirds are dissipated in inertial oscillations around the equilibrium state or 
propagated away in the form of gravity waves. This equilibrium state is an approximate 
geostrophic balance, that is, the equilibrium state is not found by solving the steady-state 
geostrophic equations. For a rotating eddy, the equilibrium is an approximate gradient 
balance between the Coriolis terms, the pressure gradient terms and the acceleration 
(curvature) terms. 

The initial disturbance can be thought of as containing a variety of length scales; the 
way the geostrophic adjustment affects each of these scales is determined by an important 
horizontal length scale, the Rossby radius of deformation (Rossby radius), given by: 


R 


d - 



(4.3) 


where c is the phase speed of the disturbance and/is the Coriolis parameter. For scales 
much smaller than the Rossby radius, rotation effects are negligible because the pressure 
gradient is much larger (short waves behave as gravity waves and disperse momentum 
rapidly from the initial disturbance), whereas for scales comparable to or larger than the 
Rossby radius, rotation effects are of the same order of magnitude of the pressure gradient 
(long waves become inertial oscillations and propagate much slower). The time scale of 
the adjustment process is the inertial period {2itlf= 19 h at 38.6°) at a distance from the 
initial disturbance comparable to the Rossby radius, but allowance has to be made for the 
oscillations about the equilibrium state. 

The considerations above can be completed with an analysis of the non- 
dimensional potential vorticity associated with the initial disturbance (Gill, 1982): 


^.(^- 5 ) =0 


(4.4) 


derived from the shallow water equations: 
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where T| is the surface elevation and H is the water column depth. 

The Rossby radius is the length scale at which, for geostrophic motion, the potential 
vorticity terms in Equation 4.4 are of the same order of magnitude; short waves are 
dominated by the vorticity term, whereas long waves are dominated by the pressure 
perturbation term ip’ = gpTj). The term including r| is related to changes in the mass 
(pressure) field; the term including ^ is related to changes in the velocity field. Therefore, 
for an initial potential vorticity that must be conserved during the adjustment, the process 
takes place differently, and the final state is different, depending on the scale of the 
disturbance compared to the Rossby radius. For scales smaller than the Rossby radius, 
where gravity is the main restoring force, the final relative vorticity dominates the potential 
vorticity, and the mass field adjusts to be in geostrophic balance with the resulting velocity 
field. For scales comparable to or larger than the Rossby radius, where rotation is 
comparable to gravity, the final mass field dominates the potential vorticity, and it is the 
velocity field that adjusts itself to a geostrophic balance with the mass field. This is 

illustrated in Figure 13. 
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Figure 13. Relation between the mass and velocity fields during the geostrophic 
adjustment, according to the length scale (L) relative to the Rossby radius of 
deformation is the initial velocity field and p is the initial mass field. 

Thus it is seen that the choice of the initial velocity field is an important issue for 
scales smaller than the Rossby radius. The density field is not the proper information with 
which to describe these scales, because it adjusts to the velocity field. Observations of the 
velocity field at small scales must be available if these scales are to be initialized in a 
circulation model. For scales comparable to or larger than the Rossby radius, description of 
the mass field is sufficient information, and the velocities adjust to reach geostrophic 
equilibrium with the observed mass field. Selecting an initial velocity field is not critical 
for these scales; the difference between starting with geostrophic initial velocities or zero 
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initial velocities is reflected only in some damping of the density field in the latter case, 
because more energy is extracted from the initial disturbance to spin-up the velocity field. 

The DPI results reveal meanders with horizontal length scales of 20 to 25 km; to 
determine whether they are influenced by the specification of the initial velocity field or 
not, an estimate of the Rossby radius is needed. In a continuously stratified fluid, with a 
vertical scale much smaller than the horizontal scale and no variations in bottom 
topography, the normal modes solution to the motion can be applied, based on a separation 
of variables between depth independent horizontal information and horizontally 
independent vertical structure. Each normal mode has its own vertical structure and 
behaves like an equivalent homogeneous system with depth Hg (equivalent depth). 

The first mode (or mode zero) is the barotropic mode with equivalent depth equal 
to the total depth of the water column. For systems predominantly baroclinic, the vertical 
structure of the higher order modes is more useful. The normal modes (barotropic and 
baroclinic) can be found by solving the following eigenvalue problem (Kundu, 1990): 


dzyfj^dz j 




(4.6) 


where is the depth dependent buoyancy frequency, \|r„ is the vertical structure function 
of any variable, is a separation constant with units of velocity (reciprocal of the 
eigenvalues) and n is the mode number. This equation was solved for the normal modes of 
vertical velocity, using an average profile for the buoyancy frequency computed from all 
the CTD casts in the cruise assuming a maximum common depth of 2000 m, with 
specification of appropriate boundary conditions for the vertical function (rigid lid at the 
surface and bottom): 

= 0, z = 0,-H. (4.7) 

The modal shapes for the first four baroclinic modes for vertical velocity are shown 
in Figure 14. The first mode contains about 20% of the total fit distributed by all the 
resolved modes (250). The other modes contribute less than 6% to the total fit. A 
comparison of Figure 14 with Figure 11 shows that the first mode describes well the region 
of maximum vertical velocities between 50 and 300 m; a combination of the remaining 
modes can account for the maximum extension of the eddy, varying from 500 to about 
1000 m. 
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Figure 14. Normalized modal shapes for vertical velocity computed from the mean 
stratification. 


Each of the resolved modes has a corresponding eigenvalue (X,=l/c„^) from which 
the associated internal Rosshy radius can be computed: 



Table 3 lists the results for the first three baroclinic modes. 


(4.8) 


Mode No. 

X 

C„ (mis) 

Rd(km) 

1 

0.235 

2.06 

23 

2 

1.070 

0.97 

11 

3 

2.474 

0.64 

7 


Table 3. Characteristics of the First Three Baroclinic Modes 

























The first internal or baroclinic Rossby radius of deformation is the horizontal length 
scale associated with coastal phenomena, such as upwelling and meandering jets, and 
eddies (Gill, 1982). These results are in good agreement with the estimates of the Rossby 
radius off the coast of Oregon, given in Kundu et al. (1975), obtained with a similar method 
to the one described above (decomposition into normal modes). 

It follows from the above that it is reasonable to assume that the vertical motion 
associated with the meanders observed in the mean flow has a vertical structure dominated 
by the first mode, which has a Rossby radius of 23 km. The length scales of the up- and 
downwelling patterns shown in Figure 9 are comparable to the Rossby radius. Thus, during 
the model initialization, it is the velocity field that primarily adjusts to the observed density 
field at the scales of the meanders. 

Therefore, the choice of the initial condition for the velocity field is not critical to 
achieve the same qualitative results from the model (meandering of the mean flow, with 
associated patches of rising and sinking motion). This is illustrated in Figure 15, which 
shows the vertical velocity at 100 m depth, obtained with zero vertical shear in the initial 
velocity field as opposed to the geostrophic shear used to obtain the results shown earlier. 
Comparing this figure with Figure 9, it is seen that the same patterns of vertical motion are 
found with the different initial velocity field, as expected. Another result that is consistent 
with the considerations above is the weaker magnitudes of the motion, which are decreased 
by more than 50 % at this level relative to the results in Figure 9. Figures 16 and 17 show 
that the vertical extent of the vertical motion cells is also somewhat smaller when the initial 
velocity field has no shear. The same qualitative agreement with the corresponding Figures 
11 and 12 is clear, but the regions of rising and sinking motion now show weaker maximum 
and minimum vertical velocities in each section. 

The question then arises as to which set of initial conditions yields results closer to 
reality, from a quantitative point of view. The results obtained when the initial velocity field 
is assumed geostrophic show a better agreement with observed currents. To set the initial 
shear velocities to zero is equivalent to saying that no a priori knowledge is available about 
the expected velocities in the eddy. However, scaling arguments can show that the eddy 
should be in a quasigeostrophic balance with the observed mass field (the Rossby number 
is of order 0.1), therefore since no detailed observations of the initial velocity field are 
available, assuming initial geostrophically balanced currents seems the better assumption. 
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Figure 16. Same as Figure 11, but obtained with an initial velocity field with no 
vertical shear. 
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3. Data related problems 


The underlying assumption in the data analysis described in section EI.B is that the 
CTD stations (plus the four SAIL stations and the satellite data) are synoptic. This is a 
common assumption in oceanography, perhaps the only viable solution to the financial 
constraints and the difficulty of data collection in the ocean. Other authors who have 
studied meso- and submesoscale phenomena (Fiekas et al., 1994; Strass, 1994; Pollard and 
Regier, 1992; Le Groupe Tourbillon, 1983) make no particular point of the fact that this 
basic assumption is made, even though the observations are collected over a number of 
days or even a week. 

The scales of the phenomena being studied are related to the more or less validity 
of the synopticity assumption. Data collected for large scale phenomena with Rossby 
numbers much smaller than unity or time scales much larger than the duration of an 
hydrographic cruise can be considered synoptic without any major argument to the 
contrary. When the scales being studied are comparable to the duration of the cruise or 
characteristic length scales of other signals, then the signal-to-noise ratio in the 
measurements becomes smaller, and aliasing and/or deformation of the scales under 
analysis may occur. 

The CTD data in this cruise took about 3 days to collect. Figure 18 shows the 
cumulative time interval between the beginning of each cast. The difference in time 
between any two stations is found by subtracting the numbers between these stations and 
the immediately preceding stations (for instance, the time interval between stations 16 and 
15 is 40.9 - 21.6 - 19.3 hours. Note that the CTD number sequence -Figure 2- does not 
correspond to the order in which the stations were occupied). The worst case of possible 
aliasing or distortion occurs when stations are close in space but separated by several hours 
in time. For example, stations 16 and 15 are apart by about 20 km in space and 19 hours in 
time. Thus an oscillation in time, uniform over the entire grid, with a period of 38 hours has 
the potential to explain all of the differences between the observed values at stations 15 and 
16. Note that the spatial aliasing problem in this case is solved by sampling theory, i.e., 
scales smaller than the separation between stations are not resolved. 

This section examines two possible sources of error in the data used in this study 
and tries to draw some, if not quantitative, at least qualitative conclusions on their effect on 
the results presented in the preceding sections. These possible sources of error are a direct 
consequence of the lack of exact synopticity in the hydrographic data. 
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Distance Alongshore (km) 


Figure 18. Cumulative time interval between CTD stations, in hours. Refer to 
Figure 2 for stations numbers (nor related to the order the stations were occupied). 

a. Internal tides 

The interaction of the barotropic tidal motions with the coastline at the shelf 
or slope causes a vertical oscillation of the pycnocline with the same range of frequencies 
of the tides. These oscillations propagate away from the region of generation as internal 
waves (internal tides); the amplitude of the internal tides depends on the strength of the 
forcing, which in turn is a function of the bottom slope, the barotropic tide—induced across 
shelf velocity and the stratification of the water colunrn (Baines, 1986). Their phase speed 
is a function of the buoyancy (AO and inertial (f) frequencies. 

Internal waves can only exist in the range of frequencies 

f<(ii<N. (4.9) 
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Hence, the main components of the barotropic tide that are forcing internal 
tides are semidiurnal. Several authors (e.g., Rosenfeld, 1990) have found that the frequency 
associated with the largest variance in internal oscillations close to the shore is the one that 
corresponds to the M 2 component (about 0.08 cycles/hour). 

Two theories try to explain the propagation of internal waves. One approach 
uses ray-like propagation, and explains important phenomena such as reflections from the 
bottom. Depending on a critical slope angle determined by the wave characteristic 
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(4.10) 


the wave will be either transmitted or reflected (Pond and Pickard, 1983). This can 
account for the concentration of energy at some depths, as waves are forced to remain at 
selected levels. Another approach decomposes the oscillation into an infinity of vertical 
dynamic modes similar to the ones depicted in Figure 14, and can explain the difference in 
the amplitudes of the oscillations at different depths. Usually the first baroclinic mode 
describes most of the variance (Rosenfeld, 1990), but other researchers have found that it 
can take a combination of several of the lowest modes to explain the vertical structure of 
the internal tide (Chiswell, 1994). 

In any case, internal wave dynamics are dependent on such a variety of 
contributions that the result is a fairly complex system that may or may not be causing 
problems in the CTD data sampling in this study. The motivation to look at this problem 
came from the finding of appreciable vertical oscillations in the upper 500 m of the water 
column, based on CTD casts taken during an anchor station about 350 km off Pt. Arena, 
thus much farther offshore than the present study region, in August 1993 (Michael Kosro, 
personal comm u nication). These oscillations had a well marked, vertically coherent 
semidiurnal character, with vertical displacements of the thermocline of about 10-15 m, 
(30 m oscillations in the isotherms at 300 m depth). Such a large and organized oscillation 
at the semidiurnal period, if present in the study region, would undoubtedly affect the 
observations in the CTD data used in this study. 

The only time series available in this study are from the current meter 
records. The instruments recorded velocity and temperature at several depths (section 
n.A.2). Unfortunately, the vertical resolution is very poor, and the period covered does not 
include the time of the CTD sampling. In the following analysis, therefore, only the last 
week of current meter data were used (spanning from 12 to 20 May 1993). Figure 19 shows 
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Figure 19. Predicted sea level at the San Francisco bay entrance, for May 1993, in 
meter. The dash-dot vertical lines mark the period of the CTD observations. 


the sea level change as predicted from tide tables for the San Francisco Bay entrance for 
May 1993; it can be seen that the barotropic M 2 forcing is roughly the same during the CTD 
cruise period and the week before (the CTD period was coincident with a new moon phase, 
thus with small M 2 amplitudes). 

A classical least-squares harmonic analysis was carried out on each of the 
18 current meter time series to evaluate the spatial phase and amplitude distribution of the 
most important tidal constituents; by doing this at all levels and comparing the phase and 
amplitude differences, an idea of the modal distribution can be obtained. The harmonic 
analysis takes the time series and fits amplitude and phase to sinusoidal basis functions with 
known frequencies. With these short time series, only a few components are resolved. 
Longer time series are not advised because of the relatively short temporal scale of the 


52 

































internal tides. In most of the cases (11 out of 18 current meters) the largest amplitude was 
associated with the M 2 . In the remaining cases, the largest amplitude was either located at 
wave-wave components, or at frequencies that are not supported by internal waves, 
without any dominant component. 

Figure 20 shows the upper 500 m of an alongshore vertical cross section of 



Figure 20. Alongshore vertical cross section of temperature ('C) from raw CTD 
data. The horizontal axis is the station number in the section. 


raw temperature data (not interpolated). In the vicinity of the eddy, the temperature 
differences observed between CTD stations are 0(10 ' ’C), being larger in the lower levels. 
The fitted amplitudes for the M 2 were in general one order of magnitude smaller. 

The attempt to fit normal modes to the M2 amplitudes was inconclusive due 
to the poor vertical resolution of the continuous records. The horizontal distribution of 
fitted phases appeared random, without any distinct pattern. The assumption of random, 
rapidly varying phase distribution of the baroclinic tides has already been used by 
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Noble et al. (1987) in their analysis of tidal currents. The results above seem to support this 
assumption in the present situation. 

The cross-correlation between temperature and horizontal velocity time 
series at each instrument was also carried out to try to identify internal wave characteristics 
at the M 2 frequency, using cross-spectra techniques described by many data analysis 
references (e.g., Bendat and Piersol, 1993). According to Figure 21, the amplitude of the 
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Figure 21. Relationship between the phase and group velocities and the velocity 
vector in a propagating internal wave (the y-axis is into the paper). C is the phase 
velocity, Cg is the group velocity, k and m are the wavenumbers (from Kundu, 
1990). 

temperature disturbance associated with internal waves is 90° out of phase with the 
u-component of velocity and in phase with the v-component. The results of the analysis of 
the coherence function and phase spectra for all the instruments are listed in Table 4. An 


Station num¬ 
ber (depth) 

Tvs. U 

Tvs. V 

Coherence 

Phase 

Coherence 

Phase 

13(100) 

0.90 

95 

0.90 

175 

9(100) 

0.80 

160 

0.80 

85 

7(100) 

0.70 

170 

0.85 

95 

8(100) 

0.90 

130 

0.90 

150 

3 (150) 

0.50 

135 

0.4 

55 

7(150) 

0.05 

70 

0.15 

45 


Table 4. Results of Currentmeter Cross-correlation Analysis for the M 2 
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Table 4. Results of Currentmeter Cross-correlation Analysis for the M 2 (Continued) 

inspection of the results shows that in only a couple of instruments are the temperature and 
velocities well correlated at the M 2 frequency, with phase values that match the 
requirements above (station 13 at 100 m and station 3 at 300 m, with an allowance made 
for phase uncertainty). A variety of situations occur for other instruments: relatively high 
values of coherence are associated with phase values that seem reversed (e.g., stations 9 and 
7 at 100 m), or relatively high values of coherence are associated with meaningless phase 
values (e.g., station 13 at 150 m, station 8 at 600 m), or temperature and velocities are 
hardly or not correlated at all at this frequency (e.g., station 7 at 150 m, station 9 at 300 m). 
Moreover, the correlation and phase values are not consistent over an appreciable depth of 
the water column. These results show a complicated vertical and horizontal pattern of 
waves at the M 2 frequency whose influence in the sampling is very likely random. 

The vertical correlation along the same water column was further 
investigated by computing the correlation and coherence functions between the 
temperature time series from the instruments in the same mooring. The values found at all 
the moorings were very small, even between instruments separated by only 50 m. The 
coherence values at the M 2 frequency were generally smaller than 0.5. The vertical 
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extension of the vertical motion patterns computed by the model and presented in section 
IV. A. 1 cannot therefore be explained by an eventual aliasing of internal tides in the CTD 
data. Internal tides may be present, as suggested by the results of the harmonic and cross¬ 
correlation analyses, but the spatial distribution of phases and amplitudes seems random 
enough to cause small contributions to the field captured in the hydrographic data. The 
chances that these randomly varying fields would become coherent and influence the 
sampling in the worst possible manner seem very small. 

Several factors may be contributing to this random distribution. The period 
of the CTD sampling was during a time of weak amplitudes in the barotropic semidiurnal 
frequency range, thus the forcing is minimal. Between 22 and 25 May there was a mix of 
relatively strong wind events and relaxations (Figure 3), causing the stratification to vary 
and affecting the maintenance of internal oscillations. Finally, the area where this study 
took place may be sufficiently close to the generation sites of the internal tides (composed 
of an infini ty of source points along an irregular shelf or slope) for the internal wave field 
to be spatially incoherent and appear random. Further offshore the wave field would 
organize itself into a vertically and horizontally coherent set of modes of propagation, 
which could explain the measurements of M. Kosro away from the coast. 

Some questions concerning the propagation of internal tides were raised by 
a study of well defined internal tide oscillations 700 km west of Santa Barbara, California 
(Levine and Richman, 1989). By performing an elliptical decomposition of the observed 
signals, they found that the predominant internal waves in that area over a period of a few 
days propagated in a direction parallel to the coast, contradicting the theories that predict 
propagation normal to the coast (Baines, 1986). They suggest that internal tides can be 
refracted by the density structure of the ocean and thus propagate in a variety of directions. 
Since the density structure in the eddy in the present study is fairly complex, these findings 
support the present assumption of small spatial coherence and highly random nature of the 
internal tide signal with a minimal influence on the observations. 

b. Eddy rotation 

This subsection is devoted to an analysis of the effect of eddy rotation on the 
solutions shown earlier. This is perhaps a more acute problem than the internal tides 
problem, because there is no way to quantify the effect of sampling a rotating body of water 
while the data is being collected. The analysis presented below assumes that there is a 
wavenumber in the eddy at least equal to two; it will be shown that a more reasonable 
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estimate is at least four, based solely on the time intervals between stations. The same 
reasoning is then applied to show that the DPI solution of a wavenumber five (Figure 7) is 
very reasonable. Translation of the eddy is neglected. According to Nof (1981), considering 
a well mixed core 300 m deep, the translation speed would be about 1.4x10'^ m/s westward, 
giving only 3 to 4 km of offshore displacement during the three days of CTD data 
collection. 

The rotation of the eddy can be estimated by taking 40 cm/s as being 
representative of the tangential horizontal speed at 30 to 40 km from the center at the 
surface, and assuming a solid body rotation for simplification. Under these premisses the 
eddy takes about 6 days, or twice the period of the CTD sampling, to complete one rotation. 
The problem here is to determine how the sub-eddy scale features have been sampled 
during the half eddy rotation. Clearly, the features themselves will not propagate at the 
same speed of the almost circular mean flow, because there exists a potential vorticity 
gradient across the mean flow that acts like an equivalent beta effect for Rossby waves, thus 
retarding the propagation speed of the features relative to the mean flow. 

A theoretical estimate of the phase speed of the sub-eddy meanders can be 
obtained by applying the results of baroclinic instability theory in the special case of a 
barotropic basic state (Holton, 1979, pp. 216), with the beta effect being replaced by its 
equivalent potential vorticity gradient, which amounts to the gradient of relative vorticity 
across the mean flow due to the flat bottom assumption. The phase speed of the 
disturbances is therefore estimated from: 


where U„ is the mean speed of the layer, k is the azimuthal wavenumber, n is the radial 
coordinate in the natural coordinate system and ^ is the azimuthally averaged relative 
vorticity given by (in natural coordinates): 


r = Y.-^ 

^ R dn 


(4.12) 


R being the distance to the center of curvature and V the tangential speed. The effect of the 
second term on the right hand side of Equation 4.11 is expected to be small for high 
wavenumbers (phase speeds of about 90 to 95% of the mean azimuthal flow), in 
accordance with the results from the numerical experiments of McCalpin (1987). 
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The mean tangential velocities were computed as a function of distance to 
the center of the eddy, by taking azimuthal averages of horizontal speed at several radial 
distances from the center, and then averaging in the vertical over the first 500 m, i.e., over 
the mean vertical extent of the sub-eddy features. The results show mean tangential speeds 
of about 6 cm/s at 10 km from the center of the eddy, increasing to about 14 cm/s at 30 km 
and 17 cm/s at 35 km, and then decreasing to about 12 cm/s at 50 km from the center. Using 
these results in Equations 4.11 and 4.12, with some radial averaging over the radial scale 
of the disturbances, gives phase speeds for the sub-eddy features of about 14 cm/s at a 
distance of 30 km from the center of the eddy, that is, about 16 days to complete one 
rotation. 

This estimate for the movement of the meanders agrees well with that 
observed during the model integrations; the animation of the horizontal velocity field at 
each time step showed an angular displacement of the sub-eddy features of about 10° about 
the center in 12 hours, or 18 days for a complete rotation. 

Therefore, during the 3 days of CTD data collection, the meanders are 
estimated to have moved about 60° cyclonically. It is clear that a wavenumber two around 
the eddy cannot create the appearance of a wavenumber five seen in the results of section 
rV.A.l by aliasing only, because the same trough (or ridge) could not be sampled three 
times during the cruise. Similarly, an analysis of Figure 18 shows that there is no time 
interval between CTD stations large enough to cause aliasing of wavenumbers smaller than 
the wavenumber five found in the DFI solution. Wavenumbers larger than five can 
theoretically be aliased, but this would imply that the length scale of the sub-eddy features 
would have to be smaller than the Rossby radius, which is the preferred length scale for the 
growth of baroclinic instabilities in a/-plane (Gill, 1982). Thus aliasing in the data due to 
meander propagation is not a cause for errors in the solution. 

Comparing the time interval between CTD stations in Figure 18 with the 
streamline and vertical motion patterns in Figure 9, and using the new information on the 
propagation of these sub-eddy features, it is seen that the maximum time interval between 
adjacent stations is 34 hours (stations 20 and 21); these two stations are located at the edges 
of a sinking cell that lays upstream of a relatively broad region of upwelling. Hence, station 
20 sampled the leading edge of the cyclonically rotating downwelling cell and station 21 
sampled the trailing edge of the same cell, 34 hours later when the cell had already rotated 
about 30°. Therefore, the net effect of the asynoptic sampling in this instance is to decrease 
the apparent width of the downwelling cell. The same reasoning applies for stations 19 and 
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22 and others. In a global sense, this explains the decrease in the horizontal scale of the 
features sampled against the sense of rotation of the eddy, as exemplified, and the increase 
in the scales of the features sampled along the sense of rotation, such as the sequences of 
stations 28 to 21. 

Other sequences of stations are not significantly affected by the propagation 
of the meanders because the time between stations was not long enough to cause an 
appreciable change in length scales; these include the sequence of stations 8, 7, 6,11, 12 
and 13 which, although sampled along the rotation, are apart by only about 13 hours total. 

It then seems reasonable to conclude that the effect of the eddy rotation, or 
meander propagation, in the CTD data is to distort the horizontal scales of some of the 
features, namely those located on the southeast side of the eddy (right side in Figure 7). The 
features on the northwest side (left side in Figure 7) must have maintained their scales with 
minimum distortion because the time interval between stations is relatively short compared 
to the angular period of the meanders. Since the troughs and ridges associated with the sub¬ 
eddy scales are about 70* apart, it would require a sequence of stations with time intervals 
on the order of several days to cause aliasing of this wavenumber. Thus, although 
wavenumber five is deformed, it is an alias-free estimate of the measured field. 

c. Interpolation parameters 

The data analysis described in section n.B makes use of the MQ-B 
interpolation method to fill the model domain with the CTD and satellite data. This method 
makes use of several interpolation parameters, the most important of which are the 
interpolation constant c and the exaggeration applied to the vertical axis (see section 
n.B.l). In addition, it was also described in section n.B that the analysis procedure used for 
the temperature and salinity fields was not the same, because the number of data points for 
temperature and salinity are different. 

At that point, the criterion for selecting the interpolation parameters was to 
achieve the best fit to the observed profiles (checked by plotting each observed profile 
against several interpolated fields) and to assure a minimum difference between extreme 
data values and extreme interpolated values. As noted, the application of these criteria 
along with the constraint in the number of observations resulted in the choice of different 
sets of parameters for temperature and salinity. 

The process of determining the sensitivity of the circulation model/DFI 
results to the set of parameters chosen requires two steps: the interpolation of the data with 
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different parameters, yielding slightly different analyzed fields each time; and running the 
model with each of the newly analyzed fields to compare results. Several such runs were 
made. It was noticed that although the analyzed fields were very sensitive to the 
interpolation constant c, the dynamically balanced velocity fields maintained the same 
qualitative agreement in the location of the meanders and vertical motion areas; the 
differences were found in the magnitude of the dynamically balanced fields (either stronger 
or weaker velocities). 

Another sensitivity test was made to estimate the influence of analyzing 
temperatures and salinities differently. In this case, the two-dimensional interpolation of 
the CTD data at each level (see section n.B.2.a) was also performed for the salinities, and 
the result was then three-dimensionally interpolated into the model domain using the same 
parameters as for the temperature interpolation. The result was a salinity field that looked 
better but had worse fits with the observed profiles and maxima and minima. The model 
DFl solution showed a weaker field of vertical motion, but in good qualitative agreement 
with the original results. 

These considerations demonstrate the more or less dependence of the results 
on the MQ-B method because of the limited data available. The choice of the parameters 
utilized in this study was dictated by criteria that guaranteed (after a long process of trial 
and error) the closest possible match between observed and analyzed fields in the regions 
where data was available. Since the computed velocities and the observed currents fi'om the 
moorings agree much better for the fields analyzed as described in II.B than they do for the 
other analysis made for sensitivity purposes (again with a limited set of data for 
comparison), it is felt that the analyzed fields from n.B and the associated results are 
correct within the given constraints. 

B. QUASIGEOSTROPHIC OMEGA EQUATION RESULTS 

The vertical motion in the study area was also diagnosed by solving the Q-vector 
form of the QG omega equation, reproduced below for convenience: 

(N^v2+y2_^)w = 2V»Q, ( 4 . 13 ) 

dz 

where 
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( 4 . 14 ) 


Q = ^ 
Po 


dx 

L3>’ 


Vp 


Vp 


with zero boundary conditions, as described in section HI.B. 

The right hand side of Equation 4.13 (known as the forcing) was computed from the 
analyzed density field, and also from the dynamically balanced density field, for 
comparison purposes. The resulting Q-vector divergence fields were smoothed with one 
pass of a smoother-desmoother filter (Kreitzberg and Perkey, 1977; Schuman, 1957), as is 
common practice in similar applications (e.g., Pauley and Nieman, 1992). An example of 
the Q-vector field in a small region of the domain at 100 m depth, obtained from the 
analyzed density field, is shown in Figure 22. The figure shows the spatially averaged 
Q-vectors, dynamic height and density fields in the area of the meander located closest to 
the offshore boundary (lower part of Figure 9; this figure can be used for exact location of 
the area depicted in Figure 22). Although a complete analysis of the dynamics of this 
meander would require information on the same fields above and below the selected depth. 
Figure 22 provides enough information to assess the behavior of the forcing of the QG 
omega equation in this area. A trough in the dynamic height field extends southward and 
offshore (roughly from the point (15,14) to the point (21,7)). The geostrophic advections 
take place parallel to the dynamic height contours. Thus upstream of the trough axis there 
is a region of relatively weak cold water advection from the center of the eddy and negative 
vorticity advection from the upstream ridge, whereas downstream of the trough there is a 
region of stronger warm water advection from outside of the eddy and positive vorticity 
advection from the trough. Cold and warm advection by the geostrophic currents are related 
to downward and upward motions, respectively (second term on the right hand side of 
Equation 3.19). Since the currents, and associated advections, generally increase toward the 
surface (the first term on the right hand side of Equation 3.19), upwelling downstream and 
downwelling upstream of the trough should be expected. Positive and negative vorticity 
advection are responsible for the propagation of the meander (in the absence of 
divergence), and their vertical structure also contributes to the vertical motion. 

The expected patterns of upwelling and downwelling cells downstream and 
upstream of the trough, respectively, are well represented by the Q-vectors, which 
generally point toward regions of rising motion. On the left of the trough the vectors are in 
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Figure 22. Field of Q-vectors computed from the analyzed density field in part of the 
domain, at 100 m. Also shown are dynamic height (solid) and density (dashed) 
contours. The X and Y coordinates are the same as in Figures 9 and 23. 

general oriented away (diverge) fi'om the region of cold water advection, where the 
dynamic height and density contours intercept at a greater angle, thus denoting the area of 
sinking motion upstream of the trough. Downstream of the trough, the Q-vectors are 
oriented with maximum convergence near the area where there is a local maximum in warm 
water advection (dynamic height and density contours intercepting almost at right angles) 
and consequent upwelling. Thus, although the effects of differential vorticity advection 
have not been considered, such an analysis is consistent with a similar analysis based only 
on the changes of relative vorticity following the motion of a parcel through the meanders 
(section IV.A.l). 

The solution for w on the left hand side of Equation 4.13 (known as the response) 
was computed with a simultaneous over-relaxation method (Haltiner and Williams, 1980). 
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Figure 23 shows the field of vertical motion at 100 m depth, obtained using the analyzed 
density field. A very good qualitative agreement with the DPI solution depicted in Figure 9 
is apparent. The two solutions reproduce the same pattern of radially aligned patches of 
sinking and rising motion, with remarkable precision as to the location of the vertical 
motions and their relative magnitudes. The patterns of rising and sinking motion deduced 
above from the Q-vectors, dynamic height and density fields (Figure 22) near the meander 
in the lower part of Figure 23 are clearly represented and in agreement with the forcing. 

Quantitatively, the vertical velocity solutions in Figures 9 and 23 are somewhat 
different. The maximum QG vertical motion is about 30% larger than the DFI solution in 
both rising and sinking areas. This difference can be due to two causes: the DFI procedure 
is appreciably affecting the initial density field, causing the resulting DFI vertical velocity 
field to be considerably different than that diagnosed from the initial density field, or/and 
the difference is due to the QG assumptions. The relative importance of the former can be 
estimated by comparing the QG omega equation results obtained from the analyzed and 
from the dynamically adjusted density fields. Figure 24 shows the latter results, which 
should be compared to those in Figure 23. It is seen that the two QG solutions are not 
significantly different, that is, the QG vertical motions diagnosed from the analyzed and 
from the DFI density fields are almost the same. This result is explained by the discussion 
of the geostrophic adjustment in section IV.A.2, where it was shown that the scales being 
initialized in the model are comparable to, or larger than, the Rossby radius. Therefore, the 
density field is hardly affected by the DFI, so the QG solutions should not be very different. 
Thus, the difference between the two solutions for w (PE and QG) is entirely due to the QG 
assumptions. 

The effects of the QG assumptions in the QG solution for w can be assessed from 
the so-called generalized omega equation, which contains the higher order effects that have 
been neglected in the QG form, namely relative vorticity, tilting, stability, ageostrophic 
advection and vertical advection effects (still considering adiabatic solutions). 
Comparisons between the generalized and the QG omega equation have been made by 
Pauley and Nieman (1992). They found that the QG omega equation solution is generally 
too large in: 1) regions of strong stability; 2) regions of large relative vorticity; 3) regions 
of sinking motion in cyclonic curvatures (rising motion in cyclonic flows can either be 
under- or overestimated). A complete analysis of the generalized omega equation is 
beyond the scope of this work, but some of the ageostrophic effects can be readily 
identified. 


63 



m 


(ajoqs ssojae) xapuf a 




Figure 23. Horizontal cross section of the quasigeostrophic vertical velocity field at 
100 m depth, in m/day, from the solution of the omega equation with the analyzed 
density field. The orientation and horizontal scales are the same as in Figure 9. Note 
that the velocity scale is different from Figure 9. 
















The stability profile used in the QG solution in Figure 23 is a spatial average of all 
the CTD casts. The maximum stability is found around 30-40 m and below that, it 
decreases sharply. The eddy has a well mixed core down to 300 m, thus stability must be 
relatively weak in the eddy at 100 m depth. 

The mean relative vorticity of the eddy is about 0.1 /. Thus neglecting relative 
vorticity versus planetary vorticity in the divergence term of the vorticity equation is by 
itself introducing a 10% error in the response of the QG omega equation. This effect tends 
to be asymmetric, causing the QG solution in general to overestimate w in regions of 
positive relative vorticity and to underestimate w in regions of negative relative vorticity. 
This effect is undoubtedly contributing to some of the overestimation of w by the QG 
solution in the meander trough shown in Figure 22. 

The QG solution for the vertical motion presented in Figure 23 seems to be strong 
in both rising and sinking motions. Thus the ageostrophic advection of vorticity and density 
that otherwise would be opposing the geostrophic advections in a cyclonic flow are missing 
in the whole domain, and not only in the areas of sinking motion. Figure 25 shows the 
ageostrophic velocity field Y^g, computed from: 




V - V 

DPI ^ g 


where V^/-/ is the horizontal velocity vector resulting from the DFI procedure, and V^is 
the initial geostrophic velocity, as well as the density and dynamic height fields for the 
same region and depth as in Figure 22. The larger vectors shown in Figure 25 represent 
ageostrophic velocities that are only about 10% of the geostrophic velocities in the same 
area. It is seen that upstream of the trough, the ageostrophic velocity vectors generally 
oppose the geostrophic flow and are associated with (weak) warm water advection and 
positive geostrophic vorticity advection, both inducing rising motion. This region is where 
a maximum of sinking motion in the QG solution is related to the cold water geostrophic 
advection. Thus, neglecting ageostrophic advection of density and geostrophic vorticity, 
which is an order of the Rossby number effect, causes the sinking motion to be 
overestimated upstream of the trough by a factor of about 10%. 

The geostrophic advection of ageostrophic vorticity upstream of the trough, which 
is also an order of the Rossby number effect, can be estimated by a careful inspection of 
Figure 25. It can be seen that the ageostrophic vectors tend to curl positively in two main 
regions, one further upstream of the trough and another just downstream of the trough. 
Between these two regions the ageostrophic vorticity is generally negative. Thus, the 
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Figure 25. Field of ageostrophic velocity vectors at 100 m depth, for the same part of 
the domain as in Figure 22. Also shown are dynamic height (solid) and density 
(dashed) contours. 

geostrophic advection of ageostrophic vorticity is positive upstream of the trough. 
Therefore, neglecting this term in the QG solution also causes an overestimate of sinking 
motion upstream of the trough by about 10%. 

Downstream of the trough, the ageostrophic velocity vectors are smaller and the 
effect of the ageostrophic advection is not as evident as it is upstream of the trough. In a 
region immediately downstream of the trough axis, the ageostrophic advection is weakly 
opposing the geostrophic warm water advection and would therefore attenuate the QG 
estimate for rising motion there. Moreover, the geostrophic flow is advecting negative 
ageostrophic vorticity into the same area, with a similar contribution for the attenuation of 
the QG rising motion estimate. Including the neglect of relative vorticity in the QG 
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divergence term, these effects all indicate an overestimate of the rising motion in and 
immediately downstream of the trough in the QG solution. Further downstream the weak 
ageostrophic advections would probably slightly enhance the rising motion computed from 
the QG omega equation, and other effects that are also order of the Rossby number (Pauley 
and Nieman, 1992) may be playing a role. This mixed effect is in accordance with the 
results of Pauley and Nieman (1992) for regions of rising motion in cyclonic curvatures. 

Another reason that could account for some of the differences between the two 
solutions is the high order numeric computations involved in the determination of the 
forcing of Equation 4.13, which makes the divergence of the Q-vectors numerically 
sensitive to small disturbances in the density and geostrophic velocity fields. 

Although the analysis above was focused on a region that shows the maximum 
discrepancy between PE and QG results, it can be generalized to the other meanders as 
well. It then seems reasonable to conclude that, in the absence of other effects, it is a 
combination of eliminating the relative vorticity from the divergence term in the vorticity 
equation, and neglecting the ageostrophic advection and vorticity, that is contributing the 
most to the quantitative differences between the QG and the DFI estimates of vertical 
velocity, with the QG values being overestimated. 

The vertical structure of the QG vertical velocity solution is shown in Figures 26 
and 27. These figures should be compared to the DFI solution in Figures 11 and 12, 
respectively. Again, a very good qualitative agreement is found between the two solutions, 
with the QG solution presenting larger values for maximum and minimum velocities, but 
with some good agreement in other areas. The vertical extent and correlation of the sinking 
and rising cells of vertical motion are also captured by the QG solution, but there is some 
variability in the deeper levels that is not found in the DFI solution. The area close to the 
northeastern boundary (right side in Figure 27) where the CTD profiles were extended in 
depth below the shallow topography also shows some unrealistic values of vertical motion. 





Figure 26. Vertical cross section of the quasigeostrophic vertical velocity field in the 
alongshore direction, in m/day, with depths in meters. The section shown is the same 
as in Figure 11 (the velocity scales are different). Note that the vertical axis is not 
linear. 
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X index (alongshore) (Min,Max) W = (-18.55,17.13) m/day 









Y index (across shore) (Min,Max) W = (-16.28,25.76) m/day 



V. CONCLUSIONS 


The main objective of this study was to achieve, for the first time, a quantitative 
description of the three-dimensional circulation in a coastal eddy, including the vertical 
velocity. The temperature and salinity data collected during one of the EEC cruises was 
thus first analyzed using multiquadric-biharmonic interpolation. The resulting density field 
was then used in a primitive equation model, to obtain three-dimensional fields of 
dynamically balanced currents and densities through a digital filter initialization (DPI) 
technique. The DPI is equivalent to a nonlinear normal mode initialization, and thus the 
dynamically balanced fields resulting from the model integrations under DPI are 
representative of the slow mode at the time of the observations. 

The results from this procedure reveal a cyclonic eddy circulation (relative vorticity 
in the order of 0.1-0.15/) about 90-100 km in diameter, with a total of five meanders each 
having a mean horizontal length scale of about 25 km, present in the flow. The maximum 
horizontal velocities in the eddy are about 40-50 cm/s at the surface, in good agreement 
with the currents observed by current meters moored in the area, and with the hull mounted 
ADCP currents estimates. 

The meanders are associated with regions of rising and sinking motion, which are 
shown to be consistent with potential vorticity considerations applied to a water parcel 
following the motion, i.e, downwelling upstream and upwelling downstream of the troughs. 
These patterns of vertical motion appear to be radially aligned, distributed on either side of 
each trough axis, and extending coherently in the vertical for about 500 m or more. The 
maximum magnitudes of the vertical velocity (either rising or sinking) occur between 50 
and 150 m depth, and are about 20-22 m/day. 

The vertical distribution of the vertical velocity amplitudes is fairly well described 
by the first dynamical mode, which has a characteristic Rossby radius of 23 km. This was 
the basis for a discussion on the process of geostrophic adjustment and on the choice of an 
appropriate initial current field. It was shown that the scales determining the main patterns 
of vertical motion are comparable to the Rossby radius, and thus it is the observed 
(analyzed) density field that determines the final state of the currents at these scales. The 
initial velocity field was therefore assumed to be geostrophic while the depth averaged flow 
was computed assuming a level of known motion (2 cm/s alongshore to the northwest at 
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2000 m depth) based on current meter data. The exact magnitude of this bottom velocity 
was found not to be significant in determining the final dynamically balanced solution. 

At this stage the eddy presents a closed circulation within 30 to 35 km from the 
center, with nearly zero average vertical velocities within this area at all levels. This means 
that the eddy is no longer acquiring relative vorticity by advection from other sources, and 
consequently it is no longer spinning up, and thus it is about to detach from the origin region 
near the coast and to drift westward (or southwestward, if the effect of the California 
Current farther offshore is taken into account). However, the effects of the poleward 
undercurrent close to the shore could not be resolved in this study because of the way the 
data were processed in the northeastern part of the domain. Therefore, allowance must be 
made for possible vorticity exchanges in the regions of the eddy that are more than 35 km 
away from the center. 

Some aspects related to the lack of exact synopticity in the hydrographic data were 
also examined. The influence of the internal tides on the sampling is thought to be minimal, 
since a classical harmonic analysis and a cross-correlation analysis of the current meter 
time series indicate that the spatial phase distribution of the internal tides was sufficiently 
random. The influence of the eddy rotation, or meander propagation, during the period 
covered by the CTD sampling was also analyzed. It was shown that a wavenumber five 
disturbance in the mean flow around the eddy could not have been aliased during sampling, 
based on the sampling pattern and on the time intervals between the hydrographic stations. 

The vertical motion was also diagnosed from the analyzed density field within the 
framework of quasigeostrophic (QG) dynamics by solving the Q-vector form of the QG 
omega equation. The results thus obtained are in close agreement with the PE/DFI results. 
The same patterns of rising and sinking motion were again found and explained using QG 
dynamics. The maximum discrepancy between the DPI and the QG vertical motion is a 
magnitude difference of about 30%, with the QG method presenting the larger values. A 
clear explanation for this was given, based on the terms neglected in the QG omega 
equation (the two most important in this study are thought to be the ageostrophic 
advections, and neglecting relative vorticity in the divergence term) and illustrated using 
an example of the actual ageostrophic velocity field in the vicinity of one of the meanders. 

From all of the results summarized above, it is concluded that a robust description 
of the circulation in the observed cyclonic eddy, including the vertical motion, has been 
achieved, complementing significantly the wealth of information that has become available 
with recent observational programs in the Californian and Oregonian coasts. The detailed 


72 


circulation in eddies such as the one studied here, i.e., with the same temporal and 
horizontal scales, and at the same stage in their life span, can be inferred from the results 
above. The methodology used in this work can also be applied to other areas of the world 
ocean where appropriate quasisynoptic data sets are available. 

Future work should be directed at the inclusion of bottom topography in the PE 
model domain. The PE results in the region of the domain closer to shore are not reliable 
because the flat bottom assumption forced some of the profiles to be extended to the 
uniform depth of 2000 m using a similarity method. The effect of a varying bottom in the 
vorticity budget will surely yield a more complete understanding of the dynamics of the 
eddy, and eliminate the undesirable effect of the extended profiles. 

Other points that can be investigated deal with the optimization of the DFI 
procedure, namely the optimization of the lowpass filter and of the filtering process, since 
it was shown that filters with fewer coefficients can better meet the constraints of the spin- 
up phase. 

Another aspect that is worth investigation is to predict the future motion of the eddy. 
Once the fields have been dynamically balanced by the DFI, then the PE model can be used 
to integrate forward and simulate the evolution of the eddy. A potential problem in this 
particular data set is the absence of information on the temporal evolution of the 
boundaries. If the boundary conditions are chosen to remain constant at inflows (as they 
were during the short initialization runs), then the eddy will not be able to drift westward 
and the general result will perhaps be an elongated eddy that remains attached to the 
boundaries where inflow takes place. Thus a careful treatment of the boundaries would be 
necessary. 
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